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I  SUMMARY 

A  three-year  program  under  AFOSR  sponsorship  has  concluded.  Four  major 
tasks  are  described:  An  atomization  experiment,  an  atomization  analysis, 
a  turbulent  reactive  flow  analysis,  and  a  vaporizing  droplet  analysis.  The 
atomization  theoretical  and  experimental  studies  consider  the  nonlinear  dis¬ 
tortion  of  a  planar  liquid  fuel  stream.  Various  modes  of  surface  instabilities 
are  identified  and  analyzed  in  detail.  Critical  values  of  parameters  related  to 
surface  tension,  stream  velocities,  disturbance  wavelength,  stream  densities, 
and  liquid  stream  thickness  are  identified  as  the  separation  between  stable 
and  unstable  interfaces.  Measurements  of  the  resulting  spray  characteristics 
are  also  presented. 

The  analysis  for  turbulent  reactive  flows  emphasizes  molecular  mixing  and 
chemical  reaction  within  a  vortical  structure.  Isolated  vortices  are  examined. 
The  vortices  are  at  the  interface  between  two  flows  of  different  composition. 
Probability  density  functions  (pdf)  and  concentration  and  velocity  profiles 
are  determined  and  compared  to  existing  experimental  data.  Temporally- 
developing  and  spatially-developing  mxing  layers  are  simulated.  Merging 
and  pairing  of  vortices  is  shown  to  be  the  major  entrainment  mechanism. 
The  larger  values  of  the  pdf  on  the  high-speed-stream  side  of  the  simulated 
mixing  layer  agree  with  experiment  but  are  shown  to  be  caused  by  a  velocity 
bias  rather  than  by  larger  mixing  rates. 

The  study  of  vaporizing  droplets  involves  Navier-Stokes  solution  for  the 
axisymmetric  flow  around  and  within  fuel  droplets.  Extensions  to  the  exist¬ 
ing  theory  have  been  made  in  two  major  areas:  (1)  variable  transport  and 
thermophysical  properties  are  considered  and  (2)  interactions  between  neigh¬ 
boring  droplets  moving  in  tandem  are  analyzed.  The  dependencies  of  drag 
coefficients,  Nusselt  numbers  and  Sherwood  numbers  upon  Reynolds  num¬ 
ber  and  transfer  numbers  are  determined.  The  vaporization  rates,  droplet 
velocities,  and  spacing  between  droplets  are  determined  as  functions  of  time. 
Critical  parameters  dividing  regions  where  droplets  collide  from  regions  where 
they  separate  are  identified. 
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II  INTRODUCTION 


This  final  report  details  the  work  performed  on  each  of  four  major  tasks  in 
the  AFOSR-sponsored  study.  A  brief  description  of  each  study  is  given  in 
this  section  with  a  more  detailed  presentation  given  in  the  next  section.  A 
reader  seriously  interested  in  any  of  the  tasks  is  strongly  encouraged  to  read 
the  publications  that  have  resulted  from  this  research.  They  are  listed  in 
Appendix  A  of  this  report. 

While  the  funding  has  ceased  and  several  publications  have  resulted  on 
each  of  the  four  tasks,  the  work  is  by  no  means  completed.  Needs  for  con¬ 
tinued  study  exist,  and  the  investigators  are  committed  to  pursuing  further 
research  in  each  of  the  major  areas. 

The  AFOSR  and,  in  particular,  Dr.  Julian  M.  Tishkoff  are  acknowledged 
for  their  support  of  these  efforts.  Their  understanding  of  the  importance  of 
these  problem  areas  in  aerospace  science  is  deeply  appreciated.  The  research 
effort  has  been  reinforced  in  that  it  relates  nicely  to  other  ongoing  research 
at  UCI  sponsored  by  other  federal  agencies  (including  another  office  in  the 
Air  Force)  and  industry. 

The  fundamental  study  has  been  productive  in  terms  of  research  results 
leading  to  a  high  number  of  publications  and  presentations,  including  several 
invited  talks  and  papers  as  indicated  in  the  appendices.  The  investigators 
were  very  active  in  spray  combustion  and  in  turbulent  flows  before  the  three- 
year  study  began;  the  position  of  UCI  in  these  fields  has,  however,  been 
strengthened  as  a  result  of  this  study,  and  we  believe  it  is  leading  nicely  to 
further  studies  on  these  problems. 

Some  major  indications  of  practical  impacts  are  already  evident.  Several 
industrial  organizations  have  sought  recently  to  align  themselves  with  UCI  in 
the  areas  of  spray  combustion  and  atomization  related  to  gas  turbine  combus¬ 
tors  and  liquid  propellant  rocket  combustion.  These  organizations  include: 
Allison  Gas  Turbine  Division  of  General  Motors,  Pratt  and  Whitney  (West 
Palm  Beach),  United  Technologies  Research  Center(UTRC),  Metrolaser,  So- 
ciete  Europeenne  Propulsion  (Vernon,  France),  and  Mitsubishi  Heavy  In¬ 
dustries  (Tokyo).  Discussions  about  a  related  spray  experiment  have  been 
conducted  with  Dr.  T.  Jackson  of  AFWRDC;  an  experiment  to  be  performed 
by  Dr.  Jackson  could  be  guided  by  theoretical  efforts  here.  Similarly,  an  ex¬ 
periment  being  conducted  by  Professor  L. Melton  of  the  University  of  Texas 
and  Dr.  M.  Winters  of  UTRC  is  being  guided  by  our  theories. 
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A  Atomization  Experiment 

It  is  well  known  that  the  efficiency  and  stability  of  many  liquid-fueled  com¬ 
bustors  can  depend  critically  upon  the  initial  sizes  and  velocities  of  the 
droplets  formed  upon  atomization  of  liquid  injected  into  the  combustor. 
In  this  program,  we  have  both  an  experimental  approach  and  a  theoret¬ 
ical/computational  approach  to  address  important  questions  about  basic 
mechanisms  for  distortion  of  the  injected  stream  and  formation  of  ligaments 
and  droplets.  Relatively  little  understanding  of  these  vital  processes  has 
existed. 

The  present  study  explores  the  basic  mechanisms  of  breakup  and  ad¬ 
dresses,  in  particular,  the  extent  to  which  properties  of  the  liquid  affect 
breakup.  Three  liquids,  selected  with  varying  values  of  viscosity  and  sur¬ 
face  tension,  are  introduced  through  a  twin-fluid,  two-dimensional  nozzle. 
A  pulsed  laser  imaging  system  is  used  to  determine  the  sheet  structure  at 
breakup,  the  distance  and  time  to  breakup,  and  the  character  of  the  ligaments 
and  droplets  formed.  Experiments  are  conducted  at  two  liquid  flow  rates 
and  five  flow  rates  of  co-flowing  air.  Two  breakup  mechanisms  are  identified. 
The  first,  termed  “cellular  breakup,”  exhibits  waves  of  equal  strength  in  the 
streamwise  and  spanwise  directions,  with  membranes  stretched  between.  In 
this  case,  the  streamwise  waves  are  broken  off  the  sheet.  The  second,  termed 
“streamwise  ligament  stretching  breakup,”  occurs  at  low  liquid  and  air  flow 
conditions.  It  is  dominated  by  streamwise  ligaments  with  thin  membranes 
of  liquid  stretched  between  them.  The  drop  size  distribution  of  the  spray  is 
characterized  by  the  breakup  mechanism.  The  cellular  breakup,  which  is  ob¬ 
tained  by  increasing  the  air  and/or  liquid  velocity,  results  in  larger  drops  at 
the  centerline  and  slightly  smaller  drops  near  the  edge  of  the  spray.  For  the 
stretched  streamwise  ligament,  the  Sauter  mean  diameter  (SMD)  at  the  cen¬ 
terline  remains  constant,  while  the  SMD  increases  significantly  away  from  the 
centerline.  Overall,  SMD  decreases  with  increasing  air/liquid  relative  veloc¬ 
ity.  The  liquid  properties  in  the  range  evaluated  do  not  change  the  breakup 
mechanisms  but  do  contribute  to  an  altered  form  of  the  breakup  in  terms  of 
the  time  and  length  scales.  The  ethanol  (lower  surface  tension  compared  to 
water  and  ethylene  glycol)  displays  smaller  scale  structures  on  the  sheet,  and 
a  greater  tendency  to  breakup  into  droplets  rather  than  form  streamwise  liga¬ 
ments.  The  ethylene  glycol  (higher  viscosity  compared  to  water  and  ethanol) 
stretches  further  before  breakup,  which  results  in  the  streamwise  ligaments 
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reaching  smaller  diameters  before  fragmenting  into  drops.  The  breakup  time 
of  the  sheet  is  dependent  on  liquid  properties  and  the  relative  velocity  of  the 
liquid  and  air.  The  ethanol  (lower  surface  tension)  demonstrates  decreased 
breakup  time,  while  the  ethylene  glycol  (higher  viscosity)  exhibits  increased 
breakup  time  relative  to  water. 

B  Atomization  Theory 

Linear  and  nonlinear  analyses  of  planar  liquid  fuel  injection  have  been  con¬ 
ducted  theoretically  and  computationally.  The  liquid-gas  interfaces  have 
been  treated  as  vortex  sheets  that  have  density,  as  well  as  velocity,  disconti¬ 
nuities  across  them.  Surface  tension  is  also  taken  into  account.  Two  general 
cases  have  been  considered:  1)  surface  disturbance  wavelength  is  very  short 
compared  to  liquid  stream  thickness  and  2)  surface  disturbance  wavelength 
is  comparable  or  much  greater  than  the  liquid  stream  thickness. 

The  method  of  vortex  dynamics  is  powerful  in  that  full  nonlinearities 
can  be  retained.  The  current  method  also  allows  for  continual  refinement  in 
order  to  obtain  the  required  resolution.  The  method  does  not  resolve  viscous 
layers,  however.  Since  the  calculations  are  two-dimensional,  distortion  of  the 
stream  but  not  break-up  is  predicted.  Nevertheless,  good  estimates  of  major 
ligament  dimensions  can  be  obtained. 

The  results  indicate  that  surface  tension  is  a  stabilizing  factor  in  the  roll¬ 
up  of  the  vortical  disturbances  at  the  liquid-gas  interfaces.  The  stability  of 
various  regions  is  identified  in  terms  of  a  dimensionless  wave  number.  For  a 
liquid  stream  of  finite  thickness,  two  types  of  waves,  shifted-symmetric  and 
symmetric  waves,  were  considered.  The  dependence  of  the  critical  dimen¬ 
sionless  wavenumber  (that  divides  stable  from  unstable  domains)  upon  the 
ratio  of  stream  thickness  to  wavelength  was  identified.  For  moderate  or  no 
surface  tension,  the  vorticity  tends  to  accumulate  near  one  point  in  the  wave, 
resulting  in  the  tendency  to  roll-up  the  interface.  Density  ratio  does  affect 
the  behavior. 

C  Turbulent  Reactive  Flows 

The  study  addresses  mixing  and  chemical  reactions  within  vortical  structures 
of  the  type  found  in  turbulent  mixing  layers.  A  series  of  investigations  has 
been  performed  in  a  sequence  that  adds  more  realism  to  the  model.  First, 
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an  isolated  vortical  structure  with  molecular  mixing  and  chemical  reaction 
in  the  windings  of  the  structure  was  analyzed.  Then  am  infinte  row  of  vor¬ 
tices  simulating  a  temporally-developing  mixing  layer  was  studied  Next,  the 
vortices  in  the  infinite  row  were  assigned  different  ages  in  order  to  simulate 
better  a  spatially- developing  mixing  layer.  Finally,  merging  of  the  vortical 
structures  was  analyzed  to  allow  for  mixing  layer  growth. 

The  major  results  indicate  a  strong  dependence  of  the  mixing  rate  on 
the  vortical  strength  as  defined  by  the  Reynolds  number  and  a  very  weak 
dependence  on  Schmidt  number.  The  probability  density  functions  (pdf’s) 
do  show  a  peak  on  the  high-speed-stream  side  of  mixing  layer  as  found  in 
experiments.  However,  this  is  not  indicative  of  a  larger  mixing  rate  on  that 
side  of  the  stream  as  interpreted  by  others;  rather,  it  indicates  that  there  is 
a  velocity  bias  in  measurement.  The  pairing  and  merging  of  vortices  is  found 
to  be  a  major  factor  in  the  entrainment  of  fluid  into  the  mixing  layer. 

D  Vaporizing  Droplet  Calculations 

The  effects  of  variable  properties  and  droplet  spacing  on  droplet  vaporization 
and  drag  are  very  important  considerations  for  practical  spray  combustion. 
This  computational  study  is  aimed  at  improving  our  basic  understanding 
and  our  computational  data  base  in  these  areas.  The  theory  has  been  ex¬ 
tended  to  consider  variable  transport  and  thermophysical  properties.  Also, 
the  behaviors  of  two  droplets  in  a  tandem  trajectory  are  considered. 

The  extension  to  variable  properties  produces  significant  quantitative  dif¬ 
ferences.  Our  study,  which  emphasizes  higher  ambient  temperatures  and 
oxidizing  environments,  has  produced  notably  different  results  for  drag  co¬ 
efficients,  Nusselt  numbers,  and  Sherwood  numbers  from  the  calculations  of 
Renksizbulut  and  co-workers  who  emphasized  moderate  ambient  tempera¬ 
tures  and  fuel  vapor  environments.  In  general,  dependencies  on  initial  droplet 
Reynolds  numbers,  ambient  temperature,  initial  droplet  temperature,  and 
fuel  characteristics  were  carefully  studied.  Correlations  for  drag  coefficients, 
Nusselt  numbers,  and  Sherwood  numbers  were  found  to  depend  primarily 
upon  Reynolds  number  and  transfer  numbers. 

Droplet  trajectories  were  analyzed  for  two  droplets  interacting  in  tandem. 
Results  indicated  that  droplet  spacing  could  increase  or  decrease  in  time,  de¬ 
pending  upon  various  factors.  The  drag  coefficient  of  the  lead  droplet  tended 
to  be  similar  to  that  of  an  isolated  droplet  until  the  droplets  were  separated  by 
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only  a  very  few  diameters.  Drag  coefficients  of  the  following  or  downstream 
droplet  were  lower  than  those  of  the  lead  droplet.  Above  some  critical  ini¬ 
tial  spacing,  the  droplet  spacing  would  increase  with  time  for  droplets  of  the 
same  initial  size  and  velocity.  The  smaller  the  downstream  droplet  relative 
to  the  lead  droplet,  the  more  likely  that  separation  of  the  droplets  would 
occur.  The  nature  of  the  heating  of  the  surface  of  the  downstream  droplet 
changes  as  the  two  droplets  become  sufficiently  close  that  the  downstream 
droplet  reaches  the  near  wake  of  the  lead  droplet. 
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Ill  RESEARCH  DISCUSSION 

A  Atomization  Experiment 

latrodiLctioa 

The  atomization  experiment  is  a  study  of  the  breakup  of  a  two-dimensional 
liquid  sheet  in  the  presence  of  co-flow  air  shear.  The  objectives  are  to  identify 
the  mechanisms  of  sheet  break-up  and  to  relate  these  mechanisms  to  the  fuel 
droplet  distribution. 

Related  Research 

Rizk  and  Lefebvre  (1980)  conducted  some  of  the  first  research  using  a  two- 
dimensional  airblast  atomizer.  Their  studies  focused  not  on  the  breakup  of 
the  sheet,  but  on  the  resulting  drop  size  distribution.  They  determined  that 
SMD  decreases  with  increasing  air  velocity,  liquid  flow  rate,  and  air/liquid 
ratio,  and  that  SMD  decreases  with  decreasing  initial  sheet  thickness.  Beck, 
Lefebvre,  and  Koblish  (1989)  developed  a  two-dimensional  airblast  atom 
izer  in  which  the  shear  air  impinges  at  a  30°  angle.  They  determined  that 
SMD  increases  with  increasing  surface  tension,  and  decreases  with  increasing 
air/liquid  mass  ratio.  Increased  viscosity  was  found  to  slightly  increase  drop 
size.  Others  have  studied  the  break-up  of  the  sheet.  Sattelmayer  and  Wittig 
(1986),  for  example,  studied  sheet  breakup  in  a  prefilming  two-dimensional 
airblast  atomizer.  They  determined  that  the  breakup  process  was  periodic, 
and  that  the  frequency  of  separation  increased  with  increasing  air  velocity. 
They  found  that  several  liquid  waves  move  to  the  edge  of  the  sheet  between 
each  event  in  which  the  liquid  is  separated.  Arai  and  Hashimoto  (1985  and 
1986)  performed  experiments  using  a  two-dimensional  nozzle  placed  inside  a 
wind  tunnel  to  provide  a  source  for  the  coflowing  shear  air.  They  noted  that 
the  instability  on  the  sheet  consisted  of  two  wave  motions,  a  stationary  wave 
in  the  streamwise  direction,  and  a  progressive  wave  in  the  spanwise  direc¬ 
tion.  They  showed  that  the  sheet  breakup  length  decreased  with  increasing 
relative  velocity,  and  increased  with  increasing  liquid  velocity  and  viscosity. 
Using  preliminary  results  with  a  two-dimensional  airblast  atomizer,  Mansour 
and  Chigier  (1989)  note  a  cellular  structure  which  forms  on  the  sheet  as  it 
disintegrates.  The  cellular  structure  appears  as  small  membranes  stretched 
between  ligaments,  and  is  attributed  to  the  different  length  scales  of  turbu¬ 
lence  in  the  shear  air  streams.  Tokuoka,  et  al.  (1985)  experimented  with  the 
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effects  of  the  angle  of  impingement  of  a  shear  air  jet  on  a  liquid  sheet.  They 
determined  that  at  higher  air  velocities,  the  sheet  breakup  length  decreased 
as  the  angle  of  impingement  increased  from  45°  to  90°. 

The  present  study  addresses  both  the  liquid  sheet  breakup,  and  the  resul¬ 
tant  droplet  size  distribution.  In  addition,  the  SMD  distribution  is  resolved 
spatially.  Parameters  of  interest  include  the  liquid  and  co-flow  air  velocities 
and  the  physical  properties  of  the  liquid. 

Experiment 


The  experimental  apparatus  consists  of  a  test  fixture  and  a  supporting 
facility.  By  way  of  summary,  the  test  fixture  is  a  two-dimensional  parallel 
flow  plane  sheet  nozzle,  a  cross  section  of  which  is  shown  in  Figure  A-l. 
The  nozzle  consists  of  a  pair  of  triangular  wedges,  which  are  separated  by 
shims  to  achieve  the  desired  thickness  for  the  liquid  sheet.  The  sheet  that  is 
produced  is  4.7  cm  wide  and  508  cm  thick.  A  Plexiglass  plenum  surrounds 
the  wedges  to  contain  the  shear  air  as  it  passes  on  either  side  of  the  sheet. 
The  air  plenum  converges  slightly  as  it  approaches  the  nozzle  tip  to  maintain 
straight  flow.  A  schematic  of  the  shear  flow  facility  is  given  in  Figure  A-2. 
The  nozzle  flow  is  directed  downwards,  and  the  liquid  is  driven  through  the 
system  by  a  gear  pump  with  a  maximum  capacity  of  157  g/sec  at  1.72  MPa. 
The  two  shear  air  lines  can  each  deliver  0.126  kg/sec. 

To  properly  image  the  atomization  of  the  liquid  sheet,  a  non-intrusive 
method  for  measuring  the  dimensions  of  the  spray  elements  was  created.  A 
pulsed  nitrogen  laser  (PRA  Model  LN1000)  was  used  to  back  light  the  sheet. 
The  short  duration  of  the  laser  pulse  (1  nsec)  is  an  attractive  feature  for 
freezing  the  sheet  dynamics.  The  UV  laser  pulse  was  expanded  and  directed 
to  impinge  on  a  sheet  paper  behind  the  liquid  sheet.  The  fluorescence  of  the 
paper  created  a  visible  back  light.  The  resultant  image  was  photographed 
using  a  35  mm  camera  (Pentax  Superprogram)  and  Professional  Ektachrome 
film  (P800/1600).  Developing  was  pushed  two  stops.  An  Aerometrics  phase 
Doppler  laser  interferometer  was  used  to  obtain  information  on  drop  size 
distribution  and  on  the  axial  velocity  of  the  spray.  The  properties  of  the 
liquids  used  for  are  given  in  Table  A-l.  Three  different  liquids  with  widely 
varying  physical  properties  (e.g.,  viscosity,  density  and  surface  tension)  were 
chosen  to  examine  the  influence  of  liquid  properties  on  sheet  breakup. 
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Table  A-l.  LIQUID  PROPERTIES 


Liquid 

Viscosity 

(Pa.sxlO) 

Surface  Tension 
(mN/m) 

Density 

(kg/m3) 

Water,  H20 

0.89 

72 

998 

Ethanol,  C2H60 

1.39 

22 

789 

Ethylene 

Glycol,  C2H602 

14.72 

48 

1114 

The  test  matrix  for  the  experiment  is  given  in  Table  A-2.  The  shear  air 
velocity  of  the  sheet  was  varied  from  0  to  60  m/sec,  with  the  same  velocity 
on  both  sides  of  the  sheet.  The  liquid  velocity  was  set  at  one  of  two  values: 
1  or  5  m/sec. 

Table  A-2.  EXPERIMENTAL  TEST  MATRIX 
Case  Velocity  (m/sec) 


Liquid 

Air  I 

Air 

11 

5 

0 

0 

12 

5 

20 

20 

13 

5 

30 

30 

14 

5 

40 

40 

15 

5 

50 

50 

16 

5 

60 

60 

21 

1 

0 

0 

22 

1 

20 

20 

23 

1 

30 

30 

24 

1 

40 

40 

25 

1 

50 

50 

26 

1 

60 

60 
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Results 

Two  mechanisms  of  sheet  breakup  have  been  identified  in  the  experi¬ 
ments.  The  first  of  these,  labelled  in  this  study  as  “cellular  breakup,”  occurs 
at  higher  relative  air-to-liquid  velocities,  but  where  the  breakup  is  well  dis¬ 
placed  from  the  nozzle  tip.  This  breakup  mechanism,  illustrated  in  Figure  A- 
3,  is  characterized  by  the  presence  of  spanwise  waves  that  are  approximately 
equal  in  strength  to  the  streamwise  vortical  waves.  As  the  sheet  is  extruded 
by  the  shearing  action  of  the  air,  the  membranes  stretch  between  the  span- 
wise  waves  and  the  streamwise  vortical  waves,  forming  cell-like  structures. 
Eventually,  the  streamwise  vortices  are  stretched  to  the  point  where  ( 1 )  the 
streamwise  waves  and  connecting  membranes  burst,  and  (2)  the  spanwise 
waves  separate  into  spanwise  ligaments.  The  small  droplets  in  the  resultant 
spray  distribution  originate  from  the  bursting  membranes.  Larger  droplets 
are  associated  with  the  breakup  of  the  streamwise  waves.  The  largest  droplets 
are  associated  with  the  fragmentation  of  the  spanwise  ligaments. 

The  second  breakup  mechanism,  referred  to  as  “stretched  streamwise  lig¬ 
ament  breakup,”  occurs  at  low  liquid  velocities.  As  shown  in  Figure  A-4 
for  Case  22,  this  mechanism  is  dominated  by  the  streamwise  vortical  waves 
in  the  breakup  process.  As  the  sheet  is  stretched  by  the  co-  flowing  air, 
the  streamwise  vortices  amplify  with  thin  liquid  membranes  stretched  be¬ 
tween.  The  sheet  is  a  corrugated  structure  with  the  membranes  stretched 
sinusoidally  about  the  plane  of  the  sheet.  A  side  view  of  this  process  is  pre¬ 
sented  in  Figure  A-4b.  As  a  result,  the  membranes  are  stretched  by  the  co- 
and  counter-rotation  of  the  bounding  vortices.  When  the  membranes  burst, 
the  liquid  film  forms  small  drops,  while  the  vortical  waves  form  streamwise 
ligaments.  These  ligaments  are  then  stretched  by  the  shearing  action  of  the 
air,  transit  to  higher  rotating  velocities,  and  finally  fragment  into  relatively 
large  drops  when  the  surface  tension  forces  can  no  longer  keep  the  ligaments 
intact. 

Measurements  of  drop  size  distribution  were  made  at  two  liquid  veloci¬ 
ties  and  three  air  velocities.  Figure  A-5a  shows  the  size  distributions  for  a 
sheet  with  a  liquid  velocity  of  5  m/sec  and  air  velocities  of  50  and  60  m/sec 
(a  40  m/sec  case  produced  drop  sizes  that  exceeded  the  capability  of  the 
instrument).  These  curves  demonstrate  the  drop  sizes  produced  as  a  result 
of  the  cellular  breakup  mechanism.  In  this  case,  the  breakup  is  uniform,  as 
the  membranes,  the  streamwise  ligaments,  and  the  spanwise  ligaments  break 
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up  along  the  centerline  of  the  spray.  The  drop  size  decreases  at  the  edges  of 
the  spray  because  the  angular  momentum  generated  by  the  cellular  breakup 
mechanism  is  not  sufficient  to  penetrate  the  shear  air  jets.  Figure  A-5b  shows 
the  size  distribution  for  a  liquid  sheet  with  a  velocity  of  1  m/sec  and  air  veloc¬ 
ities  of  40,  50,  and  60  m/sec.  The  SMD  is  constant  along  the  centerline  of  the 
spray.  Farther  from  the  centerline,  the  SMD  decreases  with  increasing  air  ve¬ 
locity.  The  minimum  that  occurs  at  the  centerline  is  more  pronounced  at  the 
lower  air  velocity  case.  This  is  a  result  of  the  stretched  streamwise  ligament 
breakup  mechanism,  in  which  the  streamwise  ligaments  form  a  corrugated 
structure,  as  described  in  Figure  A-4b.  The  thin  membranes  burst  and  are 
carried  downstream  along  the  centerline  of  the  spray,  while  the  streamwise 
ligaments  are  literally  flung  from  the  sheet  with  sufficient  angular  momentum 
to  penetrate  the  shear  air  jets  before  breakup.  Their  subsequent  fragmenta¬ 
tion  results  in  the  formation  of  the  larger  droplets  away  from  the  center  of 
the  spray. 

The  variation  in  the  liquid  properties  considered  did  not  alter  the  general 
character  of  the  two  breakup  mechanisms.  However,  the  change  in  liquid 
properties  had  a  pronounced  effect  on  the  time  and  length  scales.  Results  are 
presented  first  for  the  influence  of  liquid  properties  on  sheet  breakup  length 
and  sheet  breakup  time.  This  provides  a  global  perspective.  Secondly,  an 
examination  of  the  effect  of  liquid  properties  is  presented  for  both  the  cellular 
breakup  and  for  the  stretched  streamwise  ligament  breakup  mechanisms.  For 
the  latter,  select  conditions  are  presented  where  the  photographs  are  most 
illustrative  of  the  effects. 

Although  the  liquid  properties  of  the  three  liquids  vary,  the  liquids  were 
selected  for  their  substantial  variation  in  viscosity  and  surface  tension.  Water 
has  relatively  high  surface  tension  and  relatively  low  viscosity,  ethanol  has 
a  three- fold  lower  surface  tension  and  a  slightly  higher  viscosity  (50%)  than 
water,  and  ethylene  glycol  has  a  viscosity  much  greater  than  water  while  its 
surface  tension  is  slightly  lower  (30%).  While  it  is  likely  that  the  liquid  prop¬ 
erty  effects  presented  are  dominated  by  the  differences  in  surface  tension  and 
viscosity,  the  variation  in  other  physical  properties  and  the  limited  number 
of  fluids  considered  preclude  an  unequivocal  determination  of  the  precise  role 
of  viscosity  and  surface  tension.  As  a  consequence,  the  results  are  presented 
in  terms  of  the  generic  names  of  the  liquids  with  a  parenthetical  note  of  the 
relative  values  of  surface  tension  an.’  viscosity.  Measurements  were  made  of 
the  length  from  the  nozzle  tip  to  the  point  at  which  the  thin  film  membranes 
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fragment. 

The  results,  presented  in  Figure  A-6a,  reveal  that  the  higher  liquid  veloc¬ 
ity  leads  to  longer  sheet  breakup  lengths.  The  breakup  lengths  were  strongly 
correlated  with  the  liquid  sheet-  air  residence  time  of  interaction.  In  Figure 
A-6b,  the  data  are  replotted  in  terms  of  interaction  time  and  air-to-liquid  rel¬ 
ative  velocity,  where  interaction  time  is  defined  as  the  breakup  length  divided 
by  the  liquid  velocity  at  the  nozzle  tip.  Based  on  these  data,  it  is  apparent 
that  the  liquid  properties  have  a  notable  effect  on  the  required  air-liquid 
interaction  time  to  breakup.  The  required  interaction  time  of  the  ethanol 
is  the  shortest,  while  the  ethylene  glycol  generally  has  the  longest  required 
time.  The  former  is  the  lowest  in  surface  tension  and  density  of  the  three 
liquids,  while  the  latter  is  the  highest  in  viscosity  and  density. 

The  effect  of  liquid  properties  on  the  cellular  breakup  mechanism  is  pre¬ 
sented  in  Figure  A-7a  by  comparing  the  face  view  of  sheets  of  water  and 
ethanol  for  Case  13.  Spanwise  waves  are  formed  closer  to  the  injector  in¬ 
let  for  ethanol  (the  lower  surface  tension  liquid),  and  the  cellular  structures 
that  form  are  smaller  than  those  of  water.  It  is  also  noteworthy  that  the 
streamwise  ligaments  fragment  into  droplets  at  the  point  the  spanwise  liga¬ 
ment  breaks  off,  while  in  the  case  of  water,  the  streamwise  ligaments  persist 
well  downstream  of  the  sheet  breakup.  The  rapid  breakup  of  the  ethanol  is 
illustrated  effectively  by  a  side  view  of  sheets  of  water  and  ethanol  (Figure 
A-7b,  Case  26).  Figure  A-7c  demonstrates  the  characteristics  of  ethylene  gly¬ 
col  relative  to  the  cellular  breakup  mechanism  (Case  14).  The  high  viscosity 
of  the  ethylene  glycol  dampens  out  the  small  scale  structures  in  the  sheet 
as  it  exits  the  nozzle.  The  spanwise  waves  develop  in  a  shorter  time,  and 
appear  more  dominant  than  the  streamwise  waves.  The  cellular  structure  in 
the  ethylene  glycol  sheet  is  not  as  uniform  as  in  the  water  sheet.  After  the 
spanwise  ligaments  break  off  the  sheet,  the  streamwise  ligaments  formed  by 
the  ethylene  glycol  (higher  viscosity  fluid)  remain  intact  for  a  longer  time 
than  those  formed  by  the  water. 

The  stretched  streamwise  ligament  breakup  mechanism  is  also  affected 
by  the  liquid  properties.  Figure  A-8  is  a  face  view  of  the  three  liquids  at 
a  low  liquid  flow  condition  (Case  22).  The  membrane  formed  between  the 
streamwise  waves  is  broken  up  rapidly  into  small  drops  in  the  case  of  the  lower 
surface  tension  liquid  (ethanol).  The  streamwise  ligaments  of  the  water  break 
off  and  form  globules  and  then  break  up  further  downstream  into  drops.  The 
ethylene  glycol  stretches  further,  and  the  streamwise  ligaments  persist  longer 
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than  those  for  either  water  or  ethanol. 

Conclusions 

Over  the  range  of  conditions  and  liquid  properties  examined,  the  following 
conclusions  can  be  made: 

•  The  sheet  breakup  follows  one  of  two  basic  breakup  mechanisms:  cellu¬ 
lar  breakup  or  stretched  streamwise  ligament  breakup.  Cellular  breakup 
occurs  at  higher  liquid  and  relative  velocities.  It  is  characterized  by  the 
presence  of  spanwise  waves  that  are  approximately  equal  in  strength  to 
the  streamwise  vortical  waves,  forming  cell-like  structures.  Stretched 
streamwise  ligament  breakup  occurs  at  low  liquid  velocities  or  where 
the  breakup  occurs  relatively  close  to  the  nozzle  at  relatively  high  air 
velocities.  It  is  characterized  by  the  dominance  of  the  streamwise  vor¬ 
tical  waves  in  the  breakup  process. 

•  Changes  in  air  and  liquid  velocity  alter  the  basic  mechanisms  of  sheet 
breakup. 

•  Sheet  breakup  length  and  time  decrease  with  increasing  relative  air/liquid 
velocity. 

•  SMD  decreases  with  increasing  air/liquid  relative  velocity.  The  cellu¬ 
lar  breakup  mechanism  produces  a  uniform  breakup  which  results  in 
a  uniform  drop  size  distribution.  The  stretched  streamwise  ligament 
breakup  mechanism  results  in  a  minimum  in  the  SMD  along  the  cen¬ 
terline. 

•  Inherent  liquid  properties  (e.g.,  surface  tension,  viscosity,  and  density) 
do  not  appreciably  affect  the  basic  mechanisms  of  sheet  breakup. 

•  Inherent  liquid  properties  do  affect  the  time  and  length  scales  of  the 
sheet. 


The  two-dimensional  airblast  atomizer  that  was  developed  for  this  study 
will  be  utilized  to  conduct  additional  unsupported  research.  Higher  liquid 
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flow  rates  will  be  studied  to  determine  at  what  point  pressure  atomization 
becomes  a  factor  in  the  breakup  and  to  add  to  the  data  on  liquid  momentum 
and  relative  velocity  effects.  The  effect  of  air/liquid  momentum  ratio  and  of 
the  relative  velocity  of  the  two  shear  air  streams  will  also  be  evaluated. 

B  Atomization  Theory 

The  mechanism  of  flame  propagation  in  a  liquid-fuel  combustor  is  directly  re¬ 
lated  to  the  spatial  spray  distribution  resulting  from  the  atomization  process. 
Therefore,  an  understanding  of  this  highly  nonlinear  problem  is  a  crucial  step 
in  the  design  of  efficient  spray-combustion  systems.  Atomization  is  typically 
characterized  by  the  formation  of  droplets  with  an  average  diameter  much 
smaller  than  the  thickness/diameter  of  the  liquid  sheet/jet.  Our  approach  to 
the  study  of  this  process  is  characterized  by  considering  a  thin  but  finite  sheet 
of  one  fluid  bounded  on  each  side  by  a  surface  of  discontinuity  (vortex  sheet) 
separating  it  from  another  fluid  of  different  density  and  including  the  effect 
of  surface  tension.  When  an  interface  is  subject  to  an  initial  disturbance  of 
infinitesimal  amplitude,  well  known  results  from  linear  theory  can  be  used  to 
predict  a  dimensionless  cut-off  wave  number  which  divides  the  wave  number 
spectrum  into  a  region  of  stable  travelling  waves  and  a  region  of  exponen¬ 
tially  growing  waves  (Lamb,  1932).  The  dimensionless  cut  off  wave  number 
is  given  by: 

(1) 

where  K  is  the  wave  number,  <t  is  the  surface  tension,  Af7  is  the  velocity 
discontinuity  across  the  interface,  and  pi  and  p2  are  the  densities  of  the 
two  fluids  p  =  pi/p2  where  pi  is  the  gas-phase  density.  Exponential  growth 
corresponds  to  W  <  1  and  the  optimum  gas-phase  density.  Exponential 
growth  occurs  for  W  =  2/3.  The  dimensionless  disturbance  evolution  is 
given  by  rj  ~  exp{uit)  where: 

U  =  i«ei^£±+2«P--~WM°S  (2) 

P  2  +  Pi  P  2  +  pi 

Here  i  is  the  imaginary  unit  so  that  the  first  term  on  the  right  hand  side  gives 
a  travelling  wave  and  the  second  one  gives  the  exponential  growth.  The  linear 
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result  ceases  to  be  valid  when  the  disturbance  amplitude  is  no  longer  small 
and  later-time  evolution  of  such  a  disturbance  cannot  be  predicted. 

In  the  initial  calculations  reported,  only  one  side  of  the  liquid  sheet  was 
considered.  This  can  be  interpreted  as  the  case  in  which  the  sheet  thickness 
is  much  larger  than  the  wavelength.  The  situation  under  study  is  depicted 
in  Figure  B.l.  The  frame  of  reference  is  moving  with  a  velocity  Ure{  = 
(Ui  +  U2)/2  where  U\  and  U2  are  the  velocities  of  the  upper  and  lower  fluid, 
respectively.  The  velocity  discontinuity  is  A U  =  U\  —  U2. 

In  the  vortex-method  approach  employed  here,  the  initially  disturbed 
vortex  sheet  is  discretized  into  point  vortices  which  then  proceed  to  move 
under  the  location  of  the  interface  at  successive  time  steps.  The  location  of 
the  vortices  at  any  time  is  found  by  solving  the  equation: 


dt 


N 


—  ^0 
J=I 


(3) 


where  Uij  is  the  velocity  induced  by  vortex  j  on  vortex  i  and  the  summation 
is  taken  over  all  vortices.  The  velocity  will,  in  general,  be  a  function  of 
the  distance  between  vortices  i  and  j,  and  of  the  strength  (circulation)  of 
vortex  j. 

The  development  of  the  vortex  discretization  method  and  the  difficulties 
associated  with  it  may  be  found  in  Fink  and  Soh  (1978)  and  Van  de  Vooren 
(1980)  and  in  the  references  cited  therein.  As  a  result  of  the  density  discon¬ 
tinuity  across  the  interface  and  the  presence  of  the  surface  tension  associated 
with  it,  the  circulation  of  an  individual  vortex  does  not  remain  constant  in 
time.  An  evolution  equation  for  the  circulation  associated  with  each  vortex 
may  be  derived.  In  dimensionless  form: 


dru 

dt 


=  2  At 


1  ,  d 
+  4  AtT> 


2  8  k 

A -  We'^Asi 
os 


(4) 


where  7,  is  the  circulation  of  vortex  i,  At  =  ( p2  —  pi)/(p2  +  pi)  is  the  Atwood 
number,  We  =  (pi  +  p2)^(AU)2 /a  is  the  Weber  number,  As,  is  the  segment 
of  the  interface  associated  with  vortex  i,  k,  is  the  curvature  of  the  interface 
at  the  location  of  vortex  i,  and  s  is  the  coordinate  along  the  interface.  Here 
we  note  that  the  first  two  terms  on  the  right  hand  side  of  (4)  are  due  to  the 
density  jump  across  the  interface,  and  the  last  term  is  due  to  the  surface 
tension  a.  Realize  that  we  are  concerned  with  large  velocity  discontinuities 
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and  small  wave  lengths  A  so  that  the  Froude  number  Fr  =  (/XU)2/(Xg)  is 
large  compared  to  one,  and  gravity  effects  are  therefore  negligible. 

The  evolution  of  an  initially  sinusoidal  disturbance  is  determined  by  track¬ 
ing  the  motion  of  the  discretized  vortices.  This  is  done  by  simultaneously 
integrating  the  equations  for  the  vortex  location  and  circulation  (vorticity). 
At  t  =  0  the  interface  is  discretized  into  N  vortices  equally  spaced  along  the 
interface  coordinate  s.  As  time  evolves,  some  portions  of  the  interface  become 
stretched  while  other  portions  become  compressed.  Consequently,  vortices 
located  in  stretching  parts  of  the  interface  are  further  separated  while  those 
located  in  compressing  parts  of  the  interface  tend  to  agglomerate.  In  order  to 
handle  this  behavior,  the  technique  employed  here  consists  in  rediscretizing 
the  interface  after  each  integration  step  by  means  of  linear  interpolation  and 
introducing  a  new  vortex  whenever  the  vortex  separation  measured  along  s 
increases  beyond  an  established  percentage  of  the  initial  separation.  Once 
the  interface  is  rediscretized  by  redistributing  the  vortices,  the  individual 
vortex  strengths  are  recalculated  by  interpolation  of  the  functional  variation 
of  the  strength  versus  the  interface  coordinate. 

The  evolution  of  an  initially  small  sinusoidal  disturbance  such  as  the  one 
described  above  is  controlled  by  three  dimensionless  parameters  which  are: 
the  density  ratio,  p  (or  equivalently,  the  Atwood  number  At  —  (1  —  p)/{l  +  p)\ 
the  dimensionless  wave  number,  W  (or  equivalently,  the  Weber  number  We  ~ 
87r/(jy(l  —  A2)))  and  the  dimensionless  initial  amplitude,  e.  The  effect  of 
viscosity,  which  is  an  important  dampening  mechanism,  is  not  considered  in 
this  calculation.  In  all  the  cases  considered  in  the  present  investigation,  the 
dimensionless  initial  amplitude  is  e  =  0.025  which  falls  within  the  range  of 
validity  of  the  linear  theory.  In  fact,  linear  theory  is  used  to  obtain  the  initial 
vorticity  distribution.  We  do  not  consider  cases  in  which  the  disturbance 
amplitude  is  initially  large. 

The  effect  of  weak  surface  tension  can  be  seen  in  Figure  B.2  which  shows 
the  time  evolution  of  the  interface  shape  (column  (a))  and  the  vorticity  distri¬ 
bution  (column  (b))  for  the  case  of  p  =  1  and  W  =  0.1  .  At  a  dimensionless 
wave  number  W  —  2/3,  which  corresponds  to  the  maximum  growth  rate 
predicted  by  the  linear  theory,  the  surface  tension  is  strong  enough  to  pre¬ 
vent  roll-up.  This  occurs  because  the  high  curvature  gradient  existing  at  the 
center  of  each  wave  restricts  the  accumulation  of  vorticity  and  the  formation 
of  the  singularity.  In  Figure  B.  3  we  show  the  time  evolution  of  the  interface 
(column  (a))  and  the  vorticity  distribution  (column  s(b)  for  the  case  of  p  =  1 
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and  W  =  2/3.  As  seen  in  this  figure,  the  qualitative  features  of  the  evolution 
are  entirely  different  from  those  corresponding  to  zero  or  weak  surface  ten¬ 
sion.  Initially  and  up  to  about  t  =  0.4,  the  interface  evolves  as  predicted  by 
the  linear  theory  with  vorticity  being  accumulated  at  the  midpoint.  Beyond 
this  point,  however,  the  increasing  curvature  gradient  at  the  center,  together 
with  the  high  surface  tension  coefficient,  results  in  a  strong  restoring  force 
that  prevents  the  vorticity  from  further  accumulating  at  the  midpoint.  In¬ 
stead,  two  relative  vorticity  maxima  are  formed  on  the  sides  of  the  center 
point  where  the  curvature  gradient  is  smaller.  Under  this  condition,  the  in¬ 
terface  cannot  roll-up  but  instead  stretches  with  increasing  amplitude.  (Note 
that  the  y-scale  is  different  in  the  last  frame  of  Figure  B.3.a.)  This  situa¬ 
tion  is  particularly  important  here  because  it  implies  that  the  fastest  growing 
waves  (those  for  which  W  ~  2/3)  do  not  roll-up.  The  foregoing  discussion 
indicates  that  the  unstable  portion  of  the  spectrum  (0  <  W  <  1)  is  divided 
into  at  least  two  distinct  regions,  one  corresponding  to  weak  surface  tension 
(or  longer  waves)  in  which  roll-up  occurs  and  the  other  corresponding  to 
strong  surface  tension  (or  shorter  waves)  in  which  roll-up  does  not  occur. 
In  addition,  there  is  the  stable  region  of  the  spectrum  (very  strong  surface 
tension  or  very  short  waves)  corresponding  to  W  >  1,  as  predicted  by  the 
linear  theory.  In  order  to  illustrate  the  transition  from  the  roll-up  region  to 
the  non-roll-up  region  in  the  spectrum,  Figure  B.4  shows  the  variation  with 
time  of  the  maximum  vorticity  for  selected  values  of  the  dimensionless  wave 
number  W .  For  W  =  0,  vorticity  accumulates  at  the  center  of  the  wave 
and  grows  without  bounds  for  increasing  time  to  develop  a  vorticity  singu¬ 
larity.  When  a  moderate  surface  tension  coefficient  exists  (or  for  slightly 
shorter  waves),  vorticity  still  tends  to  accumulate  and  form  a  singularity  at 
the  center  of  the  wave  although  it  does  so  at  a  lower  rate  [W  =  0.1).  At 
a  higher  value  of  the  dimensionless  wave  number  ( W  —  0.4),  a  transition 
period  is  observed  between  t  =  0.5  and  t  =  1.0  during  which  an  inflection  in 
the  maximum  vortex  strength  curve  occurs.  During  this  transition  period, 
the  surface  tension  force  is  able  to  temporarily  slow  down  further  accumu¬ 
lation  of  vorticity.  Eventually,  enough  vorticity  has  been  concentrated  and 
the  vorticity  runaway  occurs.  At  values  of  the  dimensionless  wave  number 
above  0.155,  the  growth  is  dominated  by  the  restoring  action  of  the  surface- 
tension  force  and,  although  some  roll-up  occurs,  a  vorticity  singularity  is  not 
formed.  At  larger  times  ( t  =  2),  the  maximum  vortex  strength  is  higher  for 
W  =  1/3  than  for  W  =  2/3  because  the  former  case  involves  substantially 
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more  stretching  of  the  interface.  The  effect  of  different  density  ratios  is  also 
investigated  for  W  =  2/3,  since  this  corresponds  to  the  optimum  growth  rate. 
Figure  B.5  illustrates  the  effect  of  the  density  ratio  on  the  maximum  vortex 
strength  for  W  =  2/3.  There  is  an  optimum  density  ratio  that  results  in  the 
fastest  accumulation  of  vorticity.  This  may  be  explained  by  the  fact  that, 
at  the  optimum  dimensionless  wavelength  number  {W  =  2/3),  the  Weber 
number  increases  with  decreasing  density  ratio  thus  producing  the  strongest 
surface-tension  force  when  p\  ~  p2-  These  cases  are  surface  tension  stabi¬ 
lized  (shorter  waves).  As  p  — ►  0  the  importance  of  the  surface-tension  force 
is  diminished  (longer)  waves),  but  the  gradual  shift  of  the  center  of  vortic¬ 
ity  accumulation  results  in  stable  traveling  waves.  These  cases  are  inertially 
stable. 

The  next  effort  was  aimed  at  extending  the  mathematical  and  numerical 
model  described  above  (see  also  Rangel  and  Sirignano,  1988)  in  order  to 
simulate  the  behavior  of  a  finite-thickness  liquid  sheet  in  contact  with  a 
stream  of  air  on  each  side  of  it.  In  this  context,  we  should  note  that  our 
previous  analysis  of  a  single  interface  represents  the  limit  of  an  infinitely  thick 
sheet  or,  more  realistically,  a  case  in  which  the  sheet  thickness  is  much  larger 
than  the  disturbance  wavelength  and  the  disturbance  amplitude.  Linear 
theory  may  be  used  also  in  the  finite-sheet  case  to  predict  the  initial  growth 
of  disturbances  on  the  liquid  sheet.  Figure  B.6  is  a  schematic  of  two  limiting 
situations  that  can  be  investigated. 

We  consider  the  surface-tension-affected  Kelvin-Helmholtz  instability  of 
a  finite-thickness  sheet  of  an  inviscid  fluid  in  contact  with  two  semi-infinite 
streams  of  a  different  fluid.  Figure  B.6  illustrates  the  two  configurations 
investigated.  In  both  cases,  the  undisturbed  thickness  of  the  sheet  is  2d  and 
the  velocity  difference  across  each  interface  is  A U.  We  considered  temporal 
instabilities  so  that  periodic  boundary  conditions  are  used  at  x  =  0  and 
x  =  A.  In  our  frame  of  reference  and  in  the  absence  of  a  disturbance,  the 
sheet  moves  with  velocity  —U  while  the  semi-infinite  streams  move  with 
velocity  U.  Two  possible  disturbance  modes  are  considered  and  illustrated 
in  Figure  B.6:  shifted  symmetric  waves  (a)  and  symmetric  waves  (b).  These 
are  defined  below. 

The  linearized  equations  for  the  inviscid  shear  flow  of  a  fluid  sheet  are 
(Squire,  1953): 

Wi=0,  V2<fo  =  0,  (5) 
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where  subscript  1  refers  to  the  fluid  on  one  side  of  the  sheet  and  subscript  2 
refers  to  the  fluid  within  the  sheet.  The  fluid  on  the  other  side  of  the  sheet 
no  d  not  be  treated  explicitly  if  symmetric  or  shifted-symmetric  conditions 
are  used.  The  kinematic  and  dynamic  conditions  at  the  interface  are: 
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is  the  dimensionless  wavenumber,  a  is  the  surface-tension  coefficient,  A  is  the 
disturbance  wavelength,  A U  is  the  velocity  difference  across  the  interface, 
p i  and  p2  are  the  densities  of  the  fluids,  and  p  =  p\j p2  is  the  density  ratio. 
The  dimensionless  variables  are  the  time  t,  the  parallel  and  normal  coor¬ 
dinates  x  and  y,  the  velocity  potential  <t>,  and  the  location  of  the  interface 
tj.  The  characteristic  length  and  time  used  in  the  nondimensionalization  are 
A  and  A/AC/,  respectively.  The  previous  equations  are  solved  by  means  of 
modal  analysis  for  the  cases  of  shifted-symmetric  and  symmetric  waves,  tn 
both  cases  the  interface  initially  located  at  y  =  d  is  perturbed  according  to 
the  relation:  tj1  =  exp[i(u>t  4-  2 7rx)],  where  u>  is  the  dimensionless  complex 
frequency.  For  the  shifted-symmetric  waves  (also  referred  to  as  antisym¬ 
metric  (Squire,  1953)  or  sinusoidal  waves),  the  interface  initially  located  at 
y  =  —d  behaves  as  rj2{x,t)  =  —  rji(x  -f  \,t),  whereas  for  the  symmetric  or 
dilational  (Hagerty  and  Shea,  1955)  waves,  the  interface  at  y  =  —d  behaves 
as  r)2(x,  t)  =  —  77!(x,t).  In  addition,  the  boundary  conditions  away  from  the 
interface  require  that  the  disturbance  vanish  as  y  — >  oo. 


Shifted-symmetric  (sinusoidal)  waves: 

In  this  case,  the  solution  of  Eqs.  (5)  that  satisfies  the  kinematic  boundary 
conditions  and  the  conditions  far  from  the  sheet  yields  the  velocity  potentials 
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(10) 


as: 
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where  h  —  d/X  is  the  ratio  of  the  half  thickness  to  the  wavelength  of  the 
undisturbed  sheet. 

The  dimensionless  complex  frequency  u>  is  obtained  from  the  dynamic 
condition,  Eq.  (8)  as: 
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Exponentially  growing  waves  occur  when  the  second  term  of  this  equation 
becomes  imaginary.  Squire  (1953)  and  Hagerty  and  Shea  (1955)  limited  their 
analysis  of  Eq.  (12)  to  the  case  of  small  density  ratios  ( p  «  1)  as  corresponds 
to  a  water  sheet  in  air  at  normal  conditions.  Squire  further  limited  his 
analysis  to  long  waves  so  that  the  approximation  tanh(2irh )  =  2 nh  could  be 
used. 

The  dimensionless  thickness,  h  =  d/X,  can  be  written  as  the  product  of 
H  and  W  in  order  to  unmask  the  presence  of  the  wavenumber  in  h.  The  new 
dimensionless  thickness  H  is  given  by: 


H  = 


d{&uy  p,P2 
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In  the  general  case  of  arbitrary  density  ratio,  the  right-hand  side  of 
Eq.  (12)  becomes  imaginary,  implying  the  existence  of  unstable  waves,  when 
the  dimensionless  wavenumber  is  below  a  critical  value  (W  <  ITC).  The 
critical  dimensionless  wavenumber  is  given  by: 


Wr  = 


1  +p 


1  +  pcoth(27r//VEc) 


(14) 


Equation  (14)  defines  a  critical  dimensionless  wavenumber  below  which 
unstable  waves  exist.  In  the  limit  of  a  very  thick  sheet  ( h  >  0.43  for  practical 
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purposes),  unstable  waves  occur  for  W  <  1.  On  the  other  hand,  as  the 
thickness  decreases,  the  critical  wavenumber  decreases  below  1.  Eq.  (14) 
also  indicates  that  there  is  a  finite  thickness  H  =  Hc  below  which  all  shifted- 
symmetric  waves  are  stable.  This  critical  thickness  corresponds  to  the  value 
of  H  that  makes  Wc  =  0  in  Eq.  (14).  For  H  <  Hc ,  the  critical  wavenumber  Wc 
becomes  negative.  To  determine  Hc  we  note  that,  for  H  in  the  neighborhood 
of  Hc,  HCW  — ►  0  and  therefore  coth(2?r//W/)  — ►  l/(2irHW)  in  Eq.  (14). 
After  simplifying  one  obtains: 


Hc  = 


1  P 
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(15) 


Symmetric  (dilational)  waves: 


The  solution  of  Eqs.(5)  that  satisfies  the  kinematic  boundary  conditions 
and  the  conditions  away  from  the  shift  for  the  symmetric  waves  yields  the 
same  velocity  potential  (Eq.  10)  for  the  fluid  above  the  sheet,  while  the 
velocity  potential  for  the  fluid  in  the  sheet  becomes: 

H  =  -i(u,  -  T)«.p[iM  +  ('«) 

The  expression  for  the  dimensionless  frequency  is  obtained  again  from  the 
dynamic  condition  and  is: 


coth(27rh)  —  p 
coth(2?rft)  +  p 


±2  7T 


(*)' 


(W[coth(27r/i)  +  p]  —  (1  +  p)  coth(27rh)}1,/2 
coth(27r/i)  -f  p 


Unstable  waves  occur  for  W  <  Wc  where: 


Wc  = 


l+P 


1  +  pta,nh.(2irHWc) 
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(18) 


Again,  the  critical  dimensionless  wavenumber  approaches  1  as  the  half- 
thickness-to- wavelength  ratio  increases  above  0.43  for  practical  purposes.  In 
contrast  with  the  shifted-symmetric  waves,  the  critical  wavenumber  of  the 
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symmetric  waves  is  greater  than  zero  for  any  finite  H  and  remains  greater 
than  zero  in  the  limit  of  an  infinitely  thin  sheet,  Wc  — ►  (1  -f  p)  as  H  — ►  0. 

Our  non-hnear  investigation  is  based,  as  it  was  in  the  single-interface  case, 
on  the  vortex  discretization  method.  It  is  described  in  detail  by  Rangel  and 
Sirignano  (1988).  The  existence  of  a  second  interface  means  that  a  second 
row  of  vortices,  located  a  distance  h  below  the  first  one,  must  be  considered. 
This  two-fold  increase  in  the  number  of  vortex  elements  translates  into  ap¬ 
proximately  a  four-fold  increase  in  computational  time.  Surface-tension  and 
density  ratio  effects  are  included  in  the  form  of  an  evolution  equation  for  the 
vortex  strength  of  each  interface.  Other  details  of  the  procedure  are  basically 
the  same  ones  described  above  and  need  not  be  repeated  here. 

Some  results  are  shown  in  Figs.  B.7.  through  B.ll.  These  figures  show 
the  evolution  of  disturbances  for  different  configurations  of  interest.  Figure 
B.7  shows  a  situation  that  approaches  the  infinitely-thick-  sheet  behavior. 
Here  the  disturbance  wavelength  is  one  half  of  the  liquid  sheet  thickness 
(note  that  the  x  and  y  scales  are  different)  and  each  interface  evolves  al¬ 
most  independently  of  the  other.  Therefore,  the  roll-up  proceeds  as  in  that 
case.  When  the  liquid  sheet  thickness  is  reduced  so  that  the  disturbance 
wavelength  is  twice  the  liquid  thickness  (Figure  B.8),  the  effect  of  the  other 
interface  begins  to  become  noticeable.  In  the  limit  of  zero  thickness  the  two 
air  streams  are  stable  because  they  have  the  same  velocity.  When  surface 
tension  is  present,  the  roll-up  features  are  gradually  eliminated.  At  the  most 
unstable  case,  the  liquid  sheet  evolves  as  illlustrated  in  Figure  B.9  ,  where 
it  may  be  argued  that  two  ligaments  will  be  produced  for  each  wavelength. 
The  symmetric  case  is  illustrated  in  Figs.  B.10  and  B.ll  for  two  different 
sheet  thicknesses.  The  roll-up  features  again  disappear  as  the  thickness  is 
decreased,  and  only  one  ligament  is  formed  for  each  wavelength  of  the  distur¬ 
bance.  Thus,  the  symmetric  mode  appears  to  yield  larger  liquid  ligaments 
and  presumably  larger  droplets  than  the  shifted-symmetric  mode. 

Further  developments  are  expected  to  include  a  parametrization  of  the 
surface-tension  and  density-ratio  effects  for  the  finite-thickness  liquid  sheet 
perturbed  in  the  symmetric  and  shifted-symmetric  modes.  Extension  to  the 
case  of  the  two  air  streams  moving  with  different  velocities  should  follow. 
By  studying  the  non-linear  evolution  of  these  disturbances,  it  is  hoped  that 
some  inferences  will  be  made  about  the  size  of  the  characteristic  ligaments 
and  effects  of  the  various  different  parameters  on  this  size. 
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C  Turbulent  Reactive  Flow 


This  task  comprises  three  parts.  In  the  first  part,  molecular  mixing  and 
chemical  reation  in  an  isolated  vortex  was  investigated  using  the  approach  of 
Marble  (1985).  The  second  part  considers  a  similar  approach  for  the  study 
of  mixing  within  a  vortex  row(approximating  a  mixing  layer),  including  the 
effect  of  the  splitter  plate  in  an  approximate  manner.  In  the  third  part,  the 
fluid  dynamics  of  vortex  pairing  and  merging  were  added  and  their  effect  on 
molecular  mixing  was  investigated. 

Mixing  and  reaction  in  the  field  of  a  vortex  were  investigated  in  the  first 
part  of  this  task.  The  results  are  included  in  the  work  of  Cetegen  and  Sirig- 
nano  (1988)  and  summarized  as  follows. 

The  characteristics  of  the  probability  density  functions  (pdf’s)  for  nonre¬ 
acting  and  reacting  scalars  in  the  field  of  a  viscous  line  vortex  were  explored 
in  this  study.  The  dependence  of  the  pdfs  on  vortex  strength  and  time  at 
Schmidt  numbers  typical  of  gases  and  liquids  were  studied.  The  pdfs  are  char¬ 
acterized  by  two  probability  peaks  at  the  unmixed  concentrations  far  from 
the  vortex  center  and  a  vortex  core  region  where  finite  probabilities  of  inter¬ 
mediate  concentrations  exist.  Other  peaks  of  mixing  are  due  to  the  molecular 
diffusion  in  the  outer  regions.  As  the  vortex  strength  increases,  the  probabil¬ 
ities  in  the  mixed  core  region  become  greater.  The  increasing  vortex  strength 
also  steepens  the  local  concentration  gradients  through  stretching,  although 
the  probability  of  the  high  concentrations  gradients  decreases  with  higher 
vortex  circulation  since  the  interface  layers  within  the  vortical  flow  winding 
become  thinner  with  increasing  strain.  The  mixedness  parameter  (given  by 
a  second  moment  of  the  pdf)  shows  that  molecular  mixing  is  enhanced  by 
increasing  vortex  strength.  The  mixedness  appears  to  be  proportional  to 
Re  =  up  to  Re  =  500.  At  large  Schmidt  numbers,  the  pdf  exhibits  two 
peaks  at  unmixed  concentrations  and  a  small  peak  in  the  core  region. 

Comparison  of  the  nonreacting  and  finite  rate  chemically  reacting  scalar 
pdfs  indicates  the  development  of  an  additional  probability  peak  at  the  de¬ 
pleted  reactant  concentration  in  the  core  region  of  the  reacting  vortex.  How¬ 
ever,  the  peaks  at  the  mixed  concentrations  for  the  nonreacting  vortex  dis¬ 
appear  for  the  reacting  vortex. 

The  study  of  mixing  of  a  passive  scalar  within  a  vortex  row  originating 
from  a  splitter  plate  has  been  described  by  Miralles-Wilhelm  et  al  (1989). 
The  major  results  are  presented  here. 
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Figure  C.l  shows  the  distortion  of  the  initially  horizontal  material  line 
separating  the  two  species  at  certain  instants  after  the  vortices  are  initiated 
(Re  =  50).  Note  that  each  one  of  the  viscous  cores  is  activated  periodically 
in  time,  so  that  the  aging  process  of  the  vortical  structures  is  considered  here. 

The  instantaneous  species  concentration  contours  for  a  case  where  Re  = 
50  and  Sc  =  1  are  presented  in  Figure  C.2.  As  time  increases  and  layers  of 
both  fluids  roll-up  into  each  other,  diffusion  takes  place  through  the  interface 
separating  the  species.  The  growth  of  the  diffusion  core  in  time  is  shown  by 
the  spread  of  intermediate  concentrations  into  the  free  streams.  The  contours 
in  both  arms  of  each  vortex  far  away  from  the  center  resemble  those  for  the 
unperturbed  flow,  i.e.  -  a  sudden  change  in  the  value  of  C.  As  one  moves 
towards  the  vortical  center,  better  mixing  is  caused  by  the  intense  interaction 
between  the  vorticity  and  the  molecular  diffusion.  The  concentration  values 
in  the  central  region  are  found  to  be  close  to  the  perfect-mixed  value  of 
C  =  0.5,  as  shown  by  the  irregularly  shaped  contours  around  the  vortical 
centers.  These  characteristics  were  also  verified  for  other  Re  and  Sc  cases. 

In  Figure  C.3,  the  concentration  profiles  for  the  case  Re  =  50,  Sc  =  1 
are  shown,  as  measured  by  a  probe  located  at  different  heights  above  the 
mixing  plane  originally  dividing  the  two  species.  Notice  that  as  the  height 
y  of  the  probe  decreases  approaching  the  mixing  plane,  a  well  mixed  region 
is  encountered  (i.e.  the  vortex  core).  It  is  in  this  region  where  the  process 
occurs  more  efficiently,  as  expected.  On  the  other  hand,  when  the  height 
of  the  probe  increases,  an  increasing  tendency  towards  non-mixing  can  be 
observed,  until  no  mixing  at  all  is  reached  far  away  from  the  vortex  center. 
In  a  frame  of  reference  attached  to  the  vortical  center,  it  has  been  shown 
by  Cetegen  and  Sirignano  (1988)  that  the  concentration  profiles  exhibit  an 
antisymmetry  so  that  C(x,  j/,<)  =  1— C(  — x,  —  y,t).  That  antisymmetry  is  not 
found  here,  for  a  frame  of  reference  attached  to  a  probe  moving  through  the 
vortical  structure.  The  oscillations  in  the  concentration  profiles  are  caused  by 
the  simultaneous  rotation  and  translation  of  the  vortical  structures  as  they 
pass  through  the  measuring  probe.  Each  concentration  jump  represents  the 
probe  crossing  the  material  interface  across  which  diffusion  is  taking  place, 
which  is  tightly  wound  at  the  time  that  the  vortical  center  reaches  the  probe. 

The  probability  density  distribution  constructed  from  these  concentration 
profiles  is  shown  in  Figure  C.4  for  the  case  where  Re  =  50,  Sc  =  1.  As  noted 
before,  the  higher  peaks  correspond  to  the  unmixed  concentrations  (C  =  0 
and  C  =  1)  far  away  from  the  vortex  center,  i.e.-the  unperturbed  flow.  Finite 
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probabilities  of  intermediate  concentrations  develop  in  the  vortex  core  region 
as  a  result  of  mixing  by  molecular  diffusion  and  vortex  interaction.  The 
probability  of  concentrations  in  the  mid-range  (0.4  -  0.6),  therefore,  increases 
in  the  core  of  the  vortex,  and  the  peaks  at  y  >  0  for  C  >  0.5  and  those 
at  y  <  0  for  C  <  0.5  are  results  of  the  vortex  motion  bringing  the  high 
concentrations  from  the  lower  plane  to  the  upper  plane  and  vice/ versa. 

The  effects  of  the  two  important  dimensionless  parameters  are  shown  in 
Figures  C.5  through  C.8.  Both  an  increase  in  the  Reynolds  number  and  a 
decrease  in  the  Schmidt  number  lead  to  a  larger  well-mixed  region  around 
the  vortical  center.  In  the  first  case,  increasing  the  vortex  strength  enhances 
the  mixing  process,  since  the  rolling-up  of  the  fluid  layers  into  each  other 
is  augmented.  On  the  other  hand,  as  Sc  increases,  the  well  mixed  region 
is  diminished  because  adjacent  stratifications  within  the  vortex  have  not 
diffused  due  to  the  low  molecular  diffusivity. 

The  numerical  results  for  the  pdf  construe '  ed  with  instantaneous  data 
taken  thoughout  the  structure  exhibi*  a  characteristic  antisymmetry  whereby 
pdf(y,C)  =  pdf(—y,l  —  C),  while  the  experimental  evidence  (Dimotakis, 
1986,  1989;  Masutani  and  Bowman,  198R)  (for  pdf’s  constructed  from  time- 
varying  data  as  the  structure  flows  past  a  probe)  shows  that  the  peaks  that 
occur  on  the  high  speed  stream  side  of  the  flow  are  greater  than  those  on  the 
low  speed  stream  side.  For  this  reason,  the  relative  difference  between  the 
pdf  at  a  given  (y,  C)  pair  and  its  antisymmetrical  counterpart  (— y,  1  —  c)  is 
plotted  in  Figure  C.9.  This  plot  shows  a  new  variable: 

D  (pdf ) high.  —  speed  (pdf)ioW  —speed  ,  ..  \ 

(pdf)  low- speed 

The  numerical  calculations  and  the  shape  of  the  graph  reveal  three  im¬ 
portant  facts.  First  of  all,  R  can  be  nonzero  as  shown  at  various  y  locations, 
and,  therefore,  antisymmetry  is  not  found.  Secondly,  the  pdf  peaks  for  the 
high  speed  stream  are  not  always  greater  than  those  for  the  low  speed  stream 
(R  can  be  a  negative  quantity).  The  other  important  issue  showed  by  this 
graph  is  that  the  absolute  value  of  R  is  much  greater  when  R  is  positive. 
This  means  that  when  the  peaks  in  the  high  speed  stream  are  greater  that 
those  in  the  low  speed  stream,  the  relative  difference  between  them  is  much 
greater  (of  the  order  of  5  times)  than  that  in  the  opposite  case  (when  the  low 
speed  stream  peaks  are  higher  that  the  high  speed  stream  ones).  This  last 
result  tends  to  agree  with  the  experimental  results  presented  by  Masutami 
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and  Bowman  (1986).  The  results  also  support  the  idea  that  these  experimen¬ 
tal  results  are  a  consequence  of  the  measurement  bias,  rather  that  a  result 
of  a  real  difference  in  mixing  rate.  This  idea  is  in  disagreement  with  the 
conclusions  reported  by  Dimotakis  (1986,  1989),  which  are  that  a  difference 
exists  in  mxiag  rate  between  the  fluid  streams. 

The  high  speed  side  of  the  vortical  structure  moves  faster  on  average 
towards  the  probe  than  the  low  speed  side,  so  that  the  probe  sees  a  different 
concentration  profile  and  a  different  pdf  from  that  seen  in  a  frame  of  reference 
attached  to  the  vortical  structure.  In  the  frame  of  reference  moving  with  the 
structure,  mixing  rates  are  identical  on  the  high  speed  and  low  speed  sides,  as 
shown  in  Figure  C.10  for  the  time  evolution  of  the  pdf.  These  time-evolving 
pdf’s  are  constructed  using  instantaneous  data  over  the  vortical  structure. 

Figure  C.ll  illustrates  the  characteristics  of  a  spatially  developing  pdf. 
When  the  probe  is  located  further  downstream,  the  vortical  structures  reach 
it  at  later  instants  of  time.  Meanwhile,  fluid  is  being  mixed,  developing  the 
peaks  at  intermediate  concentrations  in  time.  In  this  way,  as  one  moves 
downstream,  the  mixed  cores  are  larger  in  size  and  the  probabilities  of  find¬ 
ing  concentrations  in  the  mid-range  are  increased.  The  same  characteristics 
regarding  the  difference  between  the  high-speed  and  the  low-speed  streams 
are  observed  in  these  results.  However,  once  the  flow  is  developed,  the  pdf 
remains  practically  invariant  as  the  probe  is  located  further  downstream. 
Insignificant  numerical  differences  are  found  in  the  pdf  computed  when  the 
probe  is  located  at  6  or  at  9  characteristic  lengths  downstream  of  the  split¬ 
ter  plate.  It  is  believed  however  that  pairing  and  merging  of  the  vortical 
structures  would  result  in  continued  development  of  the  pdf. 

A  theoretical  explanation  of  the  lack  of  antisymmetry  of  the  concentration 
field  can  be  obtained  in  terms  of  the  conservation  equation  for  the  moment 
of  the  species  concentration.  This  moment  is  an  important  characteristic  of 
the  pdf.  Analysis  of  this  equation  yields: 


1  _ dCn  ,  1  _ dC2  1  “  dC' 

2U  dx  +2V  dy  ~2D  dx 2  +  2°  dy 2  U°  dx 


-  *'c'f  - =*■£  -  »'c'f  -  D^y  -  <2> 

In  a  frame  of  reference  attached  to  the  vortical  structure,  the  first  term  on 
the  left-hand  side  of  this  equation  can  be  antisymmetric.  In  other  frames  of 
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reference  after  a  Galilean  transformation  that  term  cannot  be  antisymmetric. 
Other  terms  are  forced  to  be  also  non-antisymmetric  in  order  to  satisfy  the 
conservation  principle.  From  the  calculations  performed,  the  dominant  terms 
are  the  third  and  fourth  terms  on  the  right-hand  side,  and  the  first  two  terms 
on  the  right-hand  side  of  the  equation,  which  are  all  non-antisymmetric.  It 
is  seen  that  the  moment  C'C'  and  therefore  the  pdf  will  be  modified  by  a 
Galilean  transformation. 

Figure  C.12  shows  the  plot  of  the  mixedness  parameter  as  a  function 
of  the  height  above  the  mixing  plane  and  the  Reynolds  number.  Based  on 
the  analysis  of  these  results,  there  seems  to  be  a  weak  dependence  on  the 
molecular  diffusion  process,  represented  by  the  Schmidt  number.  This  fact  is 
in  agreement  with  the  results  presented  by  Cetegen  and  Sirignano  (1987).  On 
the  other  hand,  the  effect  of  mixing  enhancement  by  the  vortical  structures 
appears  to  be  clear;  an  increase  in  the  vortex  strength  (i.e.  the  Reynolds 
number)  increases  considerably  the  mixedness  in  the  layer.  The  mixedness 
parameter  presented  here  is  based  upon  the  pdf  determined  by  the  simulated 
experimental  probe  and  therefore  reflects  the  measurement  bias.  That  bias 
causes  a  higher  value  of  the  mixedness  parameter  on  the  high-speed  side. 
When  the  instantaneous  data  is  used  to  construct  the  pdf,  an  antisymmetrical 
situation  results,  and  the  results  for  -f  y  and  —  y  collapse  into  one  curve.  Also, 
when  an  isolated  vortex  is  considered,  the  local  mixedness  is  higher  away  from 
the  mixing  plane,  as  compared  to  that  for  the  infinite  row  of  vortices.  For  a 
value  of  y  near  the  mixing  plane,  the  local  mixedness  is  practically  identical. 
For  all  the  cases  considered,  a  dependence  on  y/Re  is  observed. 

The  mixedness  parameter  can  be  averaged  over  y,  as  in  Cetegen  and  Sirig¬ 
nano  (1987).  Figure  C.13  shows  a  comparison  between  the  isolated  vortex 
and  the  mixing  layer  approximation  for  this  case.  The  calculations  show  a 
linear  dependence  on  Re  for  the  isolated  vortex,  and  a  dependence  on  Re0  8 
for  the  infinite  row  of  vortices.  The  first  of  these  results  is  in  agreement 
with  those  reported  by  Cetegen  and  Sirignano.  The  non-linear  behavior  of 
the  overall  mixedness  parameter  /  for  the  mixing  layer  could  be  explained 
by  the  presence  of  the  surrounding  vortices,  which  increase  the  shear  as  one 
moves  away  from  y  =  0,  as  compared  to  the  case  of  an  isolated  vortex  where 
there  is  pure  rotation. 

Figure  C.14  shows  the  streamlines  for  the  infinite  row  of  vortices.  Close  to 
the  vortical  center  the  streamlines  resemble  those  for  an  isolated  vortex,  but 
for  increasing  heights  above  the  mixing  plant  the  flow  evolves  to  a  parallel 
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flow.  In  this  way,  this  increase  in  shear  results  in  a  decrease  in  the  rolling 
up  of  fluid  layers  away  from  the  mixing  plane,  and  hence  a  decrease  in  the 
mixing  process.  This  explains  the  difference  in  mixedness  between  the  cases 
of  the  isolated  vortex  and  the  row  of  vortices  for  y  ^  0. 

The  analysis  of  vortex  pairing  and  merging,  which  represents  the  third 
part  of  this  subtask,  has  been  investigated  by  means  of  vortex  dynamics. 

In  the  method  used  here,  a  continous  distribution  of  vorticity  is  discretized 
into  vortex  lines  or  vortex  cylinders  (sometimes  named  vortex  points  and 
blobs,  respectively,  because  such  is  their  appearance  in  the  two-dimensional 
cross  section).  While  a  vortex  line  is  singular  and  has  zero  cross-  sectional 
area,  a  vortex  cylinder  has  a  finite  core  and  produces  a  velocity  field  like 
that  of  a  Rankine  vortex.  The  vortex  cylinder  is  thus  sometimes  viewed  as 
an  evolution  from  the  vortex  line  if  diffusion  were  allowed  to  act  for  some 
time.  The  use  of  finite-core  vortex  cylinders  was  introduced  by  Chorin  and 
Bernard  (1972)  as  a  means  of  reducing  the  numerical  instability  inherent  to 
this  problem.  For  an  arbitrary  array  of  vortex  lines,  or  for  a  representation 
of  a  vortex  tube  as  a  sum  of  a  large  number  of  vortex  lines,  the  line  vortex 
equations  (Lamb,  1932)  are  employed  to  follow  the  evolution  of  the  structure 
defined  by  the  vortex  array.  When  vortices  with  a  finite  cross  section  are  used, 
these  equations  are  corrected  for  vortices  falling  inside  the  core  of  another 
vortex.  In  radial  coordinates,  the  velocity  field  produced  by  a  finite-core 
vortex  (Rankine  vortex)  is: 
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When  the  continuous  and  periodic  distribution  of  vorticity  represented 
by  a  vortex  sheet  of  wavelength  A  is  approximated  by  a  discrete  distribution 
of  vorticity  with  N  line  vortices  per  wavelength,  these  equations  become: 


_1_  y,  Ar,  sinh  n(y  -  yj) 

2A  cosh  n(y  —  j/;)  —  cosh  «(i  —  Xj) 


1  ^ _ AI^  sin  k(x  —  Xj ) _ 

2A  cosh  n(y  —  y3)  —  cosh  k(x  —  Xj) 
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where  AT 3  =  7j A Sj  is  the  circulation  associated  with  the  segment  of  sheet 
A sj  represented  by  vortex  j. 

Since  the  vorticity  is  simply  convected  by  the  flow  field,  the  motion  of  the 
various  vortex  elements  is  obtained  by  integrating  the  equations 


dii 

dt 


—  Ui, 


(5) 


An  initial  distribution  of  vorticity  along  any  given  line  may  become  weak¬ 
ened  or  strengthened  depending  on  whether  the  flow  dynamics  stretches  or 
compresses  the  line.  In  order  to  retain  the  resolution  implied  by  the  ini¬ 
tial  discretization,  two  approaches  have  been  adopted  depending  on  whether 
vortex  lines  or  vortex  cylinders  are  employed.  In  the  first  approach,  using 
vortex  lines,  numerical  instability  may  occur  because  of  the  singular  character 
of  the  velocity  field  in  the  vicinity  of  any  of  the  vortices.  To  prevent  the  on¬ 
set  of  this  numerical  instability,  a  rediscretization  technique  is  used  whereby 
the  nonuniformly-distributed  vortex  elements  are  replaced  by  a  new  set  of 
uniformly-  distributed  vortices  after  each  integration  step.  This  technique 
has  been  described  in  a  study  of  interface  instability  using  vortex  dynamics 
(Rangel  and  Sirignano,  1988).  Let  with  j  =  1  to  N  be  the  locations 

of  N  vortices  defining  a  certain  surface  (a  path  in  our  two-dimensional  view). 
Setting  s  =  =  0  at  the  location  of  the  first  element  (xj,  yi),  the  coordinate 

along  the  path  is  measured  as: 


sj  =  sj  - 1  +  \/(*j  -  xi-i)2  +  (j h  ~  Vi-iY  (6) 

At  the  beginning  of  the  numerical  calculation,  the  separation  between 
any  two  consecutive  vortices  is  the  same  for  the  whole  array.  After  one  time 
integration,  the  vortices  move  to  new  space  locations  ( Xj,yj )  and,  in  gen¬ 
eral,  the  inter- vortex  separations  are  no  longer  the  same  for  all  vortices.  As 
demonstrated  by  Fink  and  Soh  (1978),  the  approximation  of  a  continuous 
distribution  of  vorticity  by  a  discrete  set  of  vortices  deteriorates  when  the 
inter-vortex  spacing  does  not  remain  the  same  for  all  vortices.  In  the  redis¬ 
cretization  scheme  employed  here,  the  s  coordinates  are  calculated  with  Eq. 
(6)  and  the  total  length,  equal  to  sn+i,  is  determined.  For  a  closed  array 
of  vortices,  there  is  no  N  -f  1  vortex  and,  instead,  the  loop  ends  at  the  first 
vortex.  The  next  step  is  to  calculate  a  new  and  uniform  inter-vortex  distance 
given  by  ss+i/N  and  determine  the  new  s- coordinate  of  each  vortex  as: 
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with  sx  =  0.  To  complete  the  procedure,  the  new  vortex  locations  (xi,yi) 
and  strengths  7,  are  determined  by  interpolation  between  the  old  locations 
and  strengths,  using  Sj  and  s,  coordinates.  In  order  to  retain  the  initial 
resolution,  a  new  vortex  is  introduced  in  the  array  whenever  the  vortex  sep¬ 
aration  distance  s^=i/N  exceeds  the  initial  separation  by  a  given  factor  (5 
percent  in  our  calculations).  Since  the  overall  circulation  must  be  conserved 
according  to  Kelvin’s  theorem,  T  =  £  AI^  =  constant,  this  condition  is  used 
to  monitor  the  accuracy  of  the  calculations. 

In  the  second  approach,  using  vortex  cylinders,  no  rediscretization  is  per¬ 
formed  because  the  presence  of  a  core  prevents  the  growth  of  the  smaller, 
numerically- related  instabilities.  Since  the  vortices  are  allowed  to  move  closer 
to  or  farther  away  from  each  other,  the  initial  resolution  does  not  remain  fixed 
in  time.  In  order  to  prevent  the  resolution  from  deteriorating  in  regions  where 
vorticity  is  being  stretched  and  in  order  to  eliminate  unnecessary  calculations 
in  regions  where  the  vorticity  agglomerates,  addition  of  new  vortices  and 
combination  of  vortices  is  performed.  Whenever  the  separation  between  any 
two  consecutive  vortices  decreases  below  a  certain  value  (typically  one  half 
of  the  initial  separation),  these  two  vortices  are  combined  and  replaced  by 
an  equivalent  single  vortex.  The  location  of  the  merged  vortex  is  determined 
by  a  circulation-averaged  mean: 


Ar(-x,'  +  Art+1xt+i 

Xm_  at,  +  Art+I  (8) 

and  the  circulation  of  the  merged  structure  is  the  sum  of  the  two  contributing 
circulations: 


Arm  =  AT,  +  ATi+1  (9) 

On  the  other  hand,  when  two  consecutive  vortices  become  too  separated  (by 
moving  apart  to  twice  their  initial  separation),  a  new  vortex  is  introduced 
halfway  between  these  two  vortices: 

x,-  +  xi+i 
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The  circulation  per  unit  length  of  the  new  vortex  is  the  average  of  the 
strengths  of  the  two  neighboring  vortices: 

_  7 »  +  7i+i 
7"  -  2 

The  circulation  of  the  introduced  vortex,  ATn,  as  well  as  that  of  the  neigh¬ 
boring  vortices,  AT,  and  ATj+i,  is  calculated  from  their  respective  strengths, 
7  and  the  lengths  As  of  the  segments  they  represent. 

We  have  investigated  the  rollup,  pairing,  and  merging  of  vortical  struc¬ 
tures  formed  during  the  evolution  of  an  inviscid  vortex  sheet  subjected  to 
a  sinusoidal  disturbance  and  one  or  more  of  its  subharmonics.  In  this  case, 
the  surface  of  discontinuity  is  discretized  into  a  large  number  of  vortex  lines 
or  vortex  cylinders,  depending  on  the  method  of  solution  employed.  A  si¬ 
nusoidal  disturbance  and  any  desired  subharmonics  are  then  applied  to  the 
interface.  Corcos  and  Sherman  (1984)  investigated  this  problem  including 
viscous  effects  using  a  finite-difference  method.  Acton  (1976)  studied  the 
evolution  of  the  inviscid  vorticity  field  using  vortex  dynamics.  Since  no  re¬ 
discretization  or  insertion  of  the  new  vortices  was  performed,  the  structure  of 
the  eddy  was  lost  once  it  was  formed.  With  the  method  employed  here,  the 
initial  resolution  can  be  maintained  to  the  extent  of  the  calculation.  In  the 
following  figures,  lengths  have  been  made  dimensionless  with  the  disturbance 
wavelength,  A,  while  time  has  been  made  dimensionless  with  a  residence  time 
defined  as  A/ A U,  where  At/  is  the  magnitude  of  the  velocity  discontinuity 
at  the  undisturbed  interface.  To  illustrate  the  difference  between  the  two 
numerical  methods  described  above,  Figure  C.15  shows  the  rollup  of  a  sin¬ 
gle  disturbance  at  a  dimensionless  time  of  2.  In  both  cases  the  calculations 
are  started  with  twenty  vortical  structures  per  wavelength.  In  Figure  Cl 5a, 
line  vortices  have  been  used  with  rediscretization  after  each  time  step.  This 
method  retains  the  initial  resolution  at  all  times. 

In  Figure  C.15b,  vortices  with  a  fixed  core  of  radius  equal  to  twice  the 
spacing  between  vortices  are  used.  This  ratio  follows  the  advice  from  Leonard 
(1980).  In  this  case,  merging  and  insertion  of  vortices  is  performed  as  ex¬ 
plained  above.  Because  the  second  method  allows  the  vortices  to  get  closer 
together,  it  results  in  a  higher  resolution  of  the  rollup  region.  Note,  how¬ 
ever,  that  the  first  method  does  not  require  the  use  of  a  core  to  eliminate 
numerical  instabilities.  Note  also  that  both  methods  predict  essentially  the 
same  aplitude  growth  which  further  indicates  that  the  main  difference  is  in 
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the  resolution  of  the  core. 

Figure  C.16  illustrates  a  situation  in  which  a  sinusoidal  disturbance  and 
its  first  subharmonic,  both  of  dimensionless  amplitude  0.025,  are  simultane¬ 
ously  imposed  on  the  interface  at  t  =  0.  These  calculations  correspond  to  the 
second  method  described  above  (using  vortex  cylinders  and  no  rediscretiza¬ 
tion),  and  the  number  of  vortex  elements  n  is  indicated  at  each  time  frame. 
Figure  C.16  shows  how  the  presence  of  the  first  subharmonic  forces  the  pair¬ 
ing  of  every  two  structures  in  the  manner  described  in  the  experiments  of  Ho 
and  Huang  (1982).  No  complete  merging  occurs  because  there  is  no  diffusion 
mechanism  present  to  spread  the  vorticity  away  from  the  interface. 

Figure  C.17  shows  the  results  for  a  case  in  which  the  second  subharmonic 
is  present  in  addition  to  the  fundamental.  In  this  case,  grouping  of  every 
three  structures  occurs  and,  most  importantly,  the  transverse  thickening  of 
the  layer  occurs  more  rapidly  since  there  are  now  three  structures  rolling  up 
together.  (Note  that  the  x  and  y  scales  in  Figure  C.17  are  different  from 
those  of  Figure  C.16.) 

Figure  C.19  shows  the  amplitude  of  the  disturbance  as  a  function  of  time 
for  the  case  of  a  single  disturbance,  a  disturbance  plus  the  first  subharmonic 
(Figure  C.16),  and  a  disturbance  plus  the  second  subharmonic  (Figure  C.17). 
The  exponential  growth  rate  predicted  by  the  linear  theory  is  valid  until  a  di¬ 
mensionless  time  of  0.8,  approximately.  From  then  on,  the  single-disturbance 
grows  slowly  to  a  maximum  and  then  actually  decreases  due  to  the  strong 
shear  generated  by  the  resulting  vortex  row  and  the  absence  of  viscous  dif¬ 
fusion.  The  cases  containing  the  first  and  second  subharmonics  continue  to 
grow  nearly  linearly  until  a  dimensionless  time  of  3.  From  this  point  on, 
the  case  containing  the  first  subharmonic  reaches  a  maximum  amplitude  and 
then  decreases  also  due  to  the  absence  of  merger  in  the  presence  of  strong 
shear.  The  amplitude  of  the  case  containing  the  second  subharmonic  contin¬ 
ues  to  grow  linearly,  although  at  a  faster  rate  induced  by  the  three  structures 
that  begin  to  roll  up  together.  Eventually,  as  a  result  of  shearing  and  the 
absence  of  viscosity,  the  amplitude  of  this  case  should  also  reach  a  maxi¬ 
mum  and  start  decreasing,  although  the  calculations  have  not  been  extended 
beyond  the  stage  shown. 

The  study  of  mixing  during  the  pairing  and  merging  of  vortical  structures 
has  concentrated  on  the  prediction  of  the  diffusion  of  a  passive  scalar  in  the 
presence  of  one,  two,  three  or  more  vortical  structures.  In  addition  to  the 
study  of  the  concentration  or  mass  fraction  field,  there  has  also  been  a  study 
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of  the  development  of  the  probability  density  function  during  the  pairing  and 
merging.  We  have  studied  situations  of  temporal  growth  in  which  the  insta- 
neous  pdf  indicates  that  mixing  is  antisymmetric.  We  have  also  studied  cases 
in  which  the  spatial  development  of  the  layer  is  approximated  in  the  temporal 
growth  studies  by  starting  the  calculations  with  a  nonuniform  distribution 
of  vorticity.  This  nonuniform  distribution  of  vorticity  is  such  that  the  size  of 
the  vortices  increases  in  the  downstream  direction,  thus  giving  some  sense  of 
spatial  development.  Upstream  and  downstream  of  the  calculation  domain, 
periodic  boundary  conditions  are  still  used  and,  therefore,  the  representation 
is  only  approximate. 

Figure  C.20  shows  the  evolution  of  the  interface  when  four  vortical  struc¬ 
tures  are  present.  The  two  on  the  right  are  larger  than  the  two  on  the  left, 
thus  the  former  are  located  downstream  and  pair  first.  In  Figure  C.21,  the 
mass  fraction  contours  of  a  passive  scalar  are  shown.  The  initial  distribu¬ 
tion  is  a  step  function,  with  zero  function,  and  one  mass  fraction  below  and 
above  the  interface,  respectively.  Figure  C.22  shows  the  instantaneous  pdf 
constructed  at  each  time  step  corresponding  to  those  of  Figure  C.21.  As  seen 
in  this  graph,  the  antisymmetry  is  lost,  but  there  is  no  indication  of  preferred 
entrainment  from  either  stream.  This  argument  is  further  strengthened  by 
examining  the  path  lines  of  Figures  C23  and  C24.  The  path  lines  indicate 
where  the  fluid  particles  actually  go  and  are  different  from  the  streamlines 
which  merely  indicate  the  direction  of  the  velocity  at  each  instant.  The  path¬ 
line  plots  of  Figs.  C.23  and  C.24  indicate  that,  to  this  level  of  resolution,  the 
same  amount  of  fluid  is  entrained  from  above  and  below  the  interface  and, 
thus,  that  no  preferred  mixing  occurs. 
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D  Vaporizing  Droplet  Calculations 


There  axe  two  tasks  in  this  study.  The  first  task  is  the  detailed  analysis 
of  a  vaporizing  droplet  in  a  convective  environment.  The  second  task  is 
the  detailed  examination  of  the  interacting  effects  between  two  vaporizing 
droplets,  which  are  moving  in  tandem  with  respect  to  the  free  stream,  in  an 
intermediate- Reynolds-number  flow.  This  research  has  been  focusing  on  the 
accurate  investigation  of  the  local  transport  processes  as  well  as  the  overall 
parameters  describing  the  system  such  as  net  drag  force,  heat  and  mass 
transfer  around  the  droplets.  The  main  objective  is  to  obtain  the  relationship 
between  the  transfer  coefficients  and  some  important  fluid  parameters,  such 
as  heat  and  mass  transfer  numbers  and  Reynolds  number,  as  well  as  geometric 
parameters,  such  as  droplet  spacing  and  size  ratio. 

We  have  considered  internal  circulation  and  transient  heating  of  the  liquid 
phase,  the  forced  convection  of  gas  phase,  the  transient  deceleration  of  the 
flow  and  variable  properties.  It  is  then  necessary  to  solve  simultaneously 
the  complete  set  of  unsteady  Navier- Stokes,  energy  and  species  equations, 
combined  with  appropriate  boundary  conditions.  The  nonlinear  and  highly 
coupled  equations  make  the  analytical  solution  almost  impossible.  We  have 
to  resort  to  an  implicit  finite-difference  numerical  algorithm. 

The  governing  equations  are  represented  in  generalized  coordinates  which 
conform  to  changing  boundaries  due  to  decreasing  droplet  radii  associated 
with  liquid  component  evaporation.  A  pressure  correction  equation  is  em¬ 
ployed  to  satisfy  indirectly  the  continuity  equation  in  the  gas  phase.  The 
liquid  stream  function  and  gas  pressure  correction  equations  are  solved  by 
the  successive  over  relaxation  (SOR)  method.  The  momentum,  energy  and 
species  equations  of  the  gas  phase,  as  well  as  the  vorticity  and  energy  equa¬ 
tions  of  the  liquid  phase,  are  solved  by  an  alternating-direction-predictor- 
corrector  (ADPC)  method.  The  non-linear  gas/liquid  interface  boundary 
equations  are  treated  by  a  quasi-linearization  technique  and  solved  directly 
by  the  inversion  of  tridiagonal  block  matrices.  The  governing  equations  of 
motion,  as  well  as  the  interface  boundary  conditions,  are  solved  sequentially 
in  an  interactive  sequence  until  convergence  is  achieved  for  each  time-step  of 
the  calculation. 

It  is  noted  that  the  generation  of  the  computational  grid  for  the  multiple- 
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droplet  case  requires  more  detailed  considerations.  The  grid  lines  have  to 
conform  to  the  shapes  of  the  droplet  surfaces  to  avoid  errors  due  to  interpo¬ 
lation.  In  the  present  calculation,  the  TOMCAT  code  is  employed  to  accom¬ 
modate  the  changing  boundary  shapes  due  to  droplet  regression,  as  well  as 
the  changing  droplet  spacing  due  to  drag  differences.  See  Raju  and  Sirignano 
(1987)  or  Chiang  and  Sirignano  (1990).  The  computational  grid  is  generated 
from  the  solution  of  a  system  of  two  Poisson  equations  in  which  the  source 
terms  are  used  to  control  the  attraction  between  grid  lines.  This  quasi-linear 
elliptic  system  of  equations  is  solved  by  finite  difference  discretization  using 
the  SOR  technique. 

Since  the  reference  frame  is  fixed  to  the  liquid  droplet  (the  lead  droplet, 
in  the  two-droplet  arrangement),  the  evaluation  of  the  drag  force  and  its 
associated  velocity  correction  for  the  gas  phase  are  incorporated  in  the  iter¬ 
ative  process.  In  order  to  maintain  a  dense  grid  distribution  at  the  droplet 
interface,  the  grid  locations  have  to  be  adjusted  at  each  time  step  to  accom¬ 
modate  droplet  surface  regression.  For  the  two-droplet  case  the  adjustment 
of  the  inter-drople*  ..  Jng  results  in  an  update  of  grid  coordinates  as  well. 
Obviously  the  r>efpcs  of  the  transformation  have  to  be  updated  whenever 
the  grid  system  is  moved. 

The  overall  procedure  iterates  the  sequential  solutions  of  the  governing 
equations  and  boundary  conditions  with  grid  and  relative  velocity  adjust¬ 
ment  until  convergence  is  achieved.  The  number  of  iterations  required  for 
convergence  in  the  pressure-correction  equation,  as  well  as  within  the  over¬ 
all  iteration  loop,  is  dependent  upon  the  physical  parameters  used  in  each 
case.  For  the  cases  of  decreased  vaporization  rates  or  small  relative  velocity 
between  two  droplets,  the  requisite  numbers  are  very  much  reduced.  After 
convergence  is  reached,  the  drag  coefficients  and  average  Nusselt  and  Sher¬ 
wood  numbers  are  evaluated  at  prescribed  time  intervals. 

Our  production  runs  do  provide  some  promising  results  toward  a  better 
understanding  of  the  vaporizing  droplets’  behaviors.  We  expect  these  re¬ 
sults  to  create  an  important  benchmark  for  the  evaluation  of  numerous  other 
simplified  models  available  today.  These  comparisons  are  essential  since  the 
results  from  the  simplified  models  should  be  used  extensively  in  practical 
applications. 

The  results  are  included  in  the  works  of  Chiang,  Raju  and  Sirignano 
[1989],  Raju  and  Sirignano  [1987,  1988],  and  Chiang  and  Sirignano  [1990] 
and  summarized  as  follows. 
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SINGLE  DROPLET  CALCULATIONS 


The  detailed  analysis  of  a  vaporizing  droplet  in  a  convective  field  with 
variable  thermophysical  properties  has  been  examined  in  this  study.  A 
parameterization  has  also  been  conducted  by  changing  the  initial  droplet 
temperature,  ambient  temperature,  initial  Reynolds  number,  fuel  type,  and 
droplet  heating  model. 

The  base  case  study  is  selected  as  a  cold  n-octane  fuel  droplet  suddenly 
injected  into  the  hot  gas  stream.  The  values  of  physical  parameters  in  the 
base  case  are  given  in  Table  D-l.  Table  D-2  summarizes  the  main  parameters 
used  for  different  cases  examined  in  this  research. 

The  results  are  presented  in  the  following  five  subsections.  The  first  three 
subsections  give  a  description  of  the  local,  as  well  as  overall,  behavior  of  the 
single  droplet  for  the  baseline  case.  The  variation  of  drag  coefficient  due  to 
differe  nt  parameters  is  given  in  the  fourth  subsection.  The  subsequent  sub¬ 
section  presents  the  possible  correlations  for  the  drag  coefficient  and  Nusselt 
and  Sherwood  numbers.  The  time  scale  used  for  the  following  discussion  is 
the  gas-phase  hydrodynamic  diffusion  time. 

Results  of  the  Global  Flow  Field 

Typical  contour  plots  of  mass  fraction,  temperature,  and  vorticity,  as  well 
as  liquid-phase  streamlines  and  the  gas-phase  velocities,  at  time  =  25  are 
presented  in  Figure  D-  l(a-e).  The  convective  effect  is  apparent  by  the  fore- 
aft  asymmetry  in  each  plot.  The  boundary  layers  surrounding  the  droplet  can 
be  easily  observed  from  the  gas-phase  velocity  distribution.  At  early  times, 
there  is  no  indication  of  flow  separation  because  the  separation  is  suppressed 
by  the  surface  movement  (which  will  be  explained  later).  At  later  time  the 
boundary  layers  only  cover  up  to  138°  in  azimuthal  position,  and  a  near 
wake  region  appears  behind  the  droplet.  The  magnitude  of  velocity  has  been 
reduced,  owing  to  the  retardation  by  the  drag  force.  The  temperature  and 
mass  fraction  contours  do  show  the  presence  of  the  boundary  layer  and  wake. 
The  diffusing  fuel-vapor  is  convected  downstream,  forming  a  concentration 
wake  at  the  rear. 

The  liquid-phase  streamlines  show  a  single  large  vortex,  roughly  resem¬ 
bling  a  Hill’s  spherical  vortex,  in  the  interior  of  droplet.  The  strength  of  cir- 
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culation  shows  a  decrease  with  time.  The  vorticity  contours  are  concentrated 
at  the  front  portion  of  the  droplet  and  then  are  convected  downstream.  At 
later  time  the  vorticity  has  been  diffused  outward  and  the  intensity  is  weak¬ 
ened,  mainly  due  to  the  reduction  in  Reynolds  number.  A  small  loop,  which 
represents  a  recirculating  eddy  attached  at  the  surface,  is  observed.  When 
flow  separation  occurs,  the  distribution  becomes  highly  asymmetric.  The 
surface  blowing  inhibits  the  appearance  of  a  small  secondary  internal  vortex 
of  opposite  direction  which  might  occur  near  the  rear  stagnation  point. 

The  detailed  droplet  heating  history  may  be  seen  from  Figure  D-  2(a-c). 
During  the  very  early  time  (<  1)  a  large  temperature  gradient  exists  near 
the  surface,  and  the  motion  of  the  liquid  phase  is  negligible.  The  energy 
transfer  mode  within  the  droplet  is  dominated  by  the  conduction.  As  the 
strength  of  internal  circulation  increases,  the  transport  mechanism  switches 
to  convection  gradually.  The  isotherms  at  time  =  6  show  the  large  tempera¬ 
ture  variation  from  the  surface  to  the  vortex  center.  The  similarity  between 
the  isotherms  and  streamlines  in  the  liquid  phase  shows  the  influence  of  in¬ 
ternal  circulation.  The  assumption  of  high  Peclet  number  used  in  a  vortex 
model  (Prakash  and  Sirignano  [1980],  Lara-Urbaneja  and  Sirignano  [1981], 
Tong  and  Sirignano  [1982])  is  reliable  at  this  stage.  When  the  gas  phase 
shear  stresses  decrease,  the  circulation  strength  lessens  gradually  such  that 
the  convection  region  becomes  smaller  and  the  region  for  conduction  grows 
again.  This  behavior  is  more  evident  for  the  lower  transfer  number  case 
where  the  Peclet  number  is  smaller. 

Local  Properties  Along  the  Droplet  Surface 

Since  the  overall  drag  coefficient,  Nusselt  and  Sherwood  numbers  of  the 
droplet  are  very  sensitive  to  the  transporting  processes  occurring  at  the 
gas/liquid  interface,  it  is  necessary  to  have  a  detailed  understanding  of  some 
important  properties  at  droplet  surface.  The  surface  shear  stress  distribu¬ 
tion  at  different  times  shown  in  Figure  D-  3(a)  indicates  that  the  shear  stress 
decreases  with  time.  This  can  be  realized  by  the  following  facts:  (1)  with 
the  onset  of  surface  motion,  the  velocity  gradient  at  the  droplet  surface  de¬ 
creases;  (2)  the  surface  blowing  effect,  which  increases  the  thickness  of  the 
boundary  layer  and  reduces  the  velocity  gradient,  is  growing  as  vaporization 
becomes  stronger;  and  (3)  as  droplet  heating  continues,  the  mass  fraction  at 
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the  surface  also  keeps  increasing,  thus  yielding  lower  values  of  viscosity  of 
the  mixture. 

The  decrease  of  strength  of  the  negative  shear  stress  is  due  to  the  decrease 
of  volume  and  strength  of  the  wake  recirculation  which  is  mainly  caused  by 
the  reduction  of  the  Reynolds  number. 

The  time  variation  of  the  surface  tangential  velocity  distribution  shown 
in  Figure  D-  3(b)  seems  to  vary  with  sinO  which  is  in  qualitative  agreement 
with  the  classical  vortex  solution.  The  surface  tangential  velocity,  originally 
zero,  is  brought  to  a  maximum  by  the  shear  stresses  and  then  decreases  as 
shear  stresses  diminish.  It  is  also  noted  that  the  angular  location  of  the 
maximum  velocity  remains  unchanged  except  during  the  final  period  when 
the  location  of  maximum  shifts  forward  of  the  6  =  90°  plane. 

The  distribution  of  vorticity  along  the  droplet  surface  during  the  transient 
development  is  very  similar  to  that  of  the  shear  stress.  The  surface  vorticity 
is  obtained  from 

v„  lav,...  „av,.e, 

=  7 3 - i - Ha- —  — I* 

fXj  on  a  a  ua  on 

The  dominating  term  is  ^Vgfi  which  is  also  the  main  term  contributing 
to  the  shear  stress.  It  is  clear  that  oj  ~  |  ~  ReJ^2.  Hence  it  is  expected  that 
the  effects  of  internal  circulation,  boundary  layer  blowing  and  reduction  of 
gas-phase  Reynolds  number  will  cause  the  surface  vorticity  to  diminish.  If 
the  flow  separation  point  is  determined  by  the  point  where  vorticity  changes 
sign,  the  results  indicate  that  separation  occurs  at  a  very  early  time.  When 
the  surface  starts  moving,  the  separation  is  suppressed.  At  later  time  the 
strength  of  surface  velocity  is  reduced,  and  the  normal  surface  blowing  ve¬ 
locity  increases.  As  a  result,  the  separation  point  moves  in  the  upstream 
direction.  This  behavior  is  different  from  that  of  the  low  transfer  number 
case  where  the  surface  velocity  is  usually  small.  The  surface  motion  can 
only  delay  the  onset  of  flow  separation.  The  flow  separation  is  observed 
throughout  the  droplet  lifetime.  The  theoretical  separation  point  for  the  liq¬ 
uid  sphere  (Clift,  Grace,  and  Weber  [1978])  moves  rearward  as  the  Reynolds 
number  decreases  with  time.  The  prediction  of  the  separation  point  for  a 
vaporizing  droplet  is  hence  totally  different  from  that  for  a  nonvaporizing 
liquid  sphere. 

Figure  D-  3(c)  shows  the  pressure  distribution  on  the  droplet  surface. 
The  recirculation  wake  dissipates  part  of  the  kinetic  energy.  As  a  result, 
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the  pressure  cannot  recover  to  the  stagnation  value.  The  pressure  profile  is 
essentially  determined  from  the  balance  between  the  diffusive  and  convective 
transport  of  vorticity.  Vaporization  enhances  the  outward  convection  of  vor- 
ticity  and  thus  causes  more  pressure  loss.  However,  the  additional  kinetic 
energy  resulting  from  surface  blowing  is  absorbed  by  the  recirculation  wake. 
This  can  result  in  some  pressure  recovery. 

The  angular  variation  of  local  Nusselt  number  at  different  times  is  shown 
in  Figure  D-  3(d).  The  Sherwood  number  behaves  similarly.  The  temperature 
(or  mass  fraction)  gradient  at  the  surface,  and  thus  Nuiocai  (or  Shiocai )  ,  is 
greatest  at  the  front  stagnation  point  and  decreases  with  polar  angle;  see  also 
Figure  D-  1(b)  and  D-  1(c).  The  minimum  does  not  necessarily  occur  aft  of 
the  separation  point,  as  for  the  case  of  low  transfer  number.  The  increase 
at  the  rear  of  the  droplet  is  caused  by  the  action  of  the  recirculating  wake. 
The  reduction  in  the  Nusselt  and  Sherwood  numbers  with  time  is  mostly 
attributed  to  the  decay  in  Reynolds  number. 

Figure  D-  3(e)  shows  the  time  variation  of  the  surface  temperature  dis¬ 
tribution.  The  surface  temperature  is  nonuniform  during  most  of  the  droplet 
lifetime.  The  rise  of  temperature  at  the  rear  portion  is  predominantly  caused 
by  the  hot  far-stream  fluid  returned  by  the  recirculating  wake.  After  droplet 
heating  begins  to  diminish,  the  surface  temperature  shows  some  degree  of 
uniformity.  However,  the  constant  surface  temperature  assumption  is  im¬ 
proper  when  the  droplet  heating  persists. 

The  transient  development  of  surface  normal  velocity  is  shown  in  Figure 
D-  3(f).  The  gas  surface  normal  velocity  is  primarily  determined  by  the 
rate  of  evaporation.  Hence,  the  normal  velocity  increases  when  the  surface 
temperature  increases  with  time.  As  the  surface  temperature  approaches  the 
wet-bulb  temperature,  the  fuel-vapor  mass-fraction  at  the  surface  remains 
constant  and  the  boundary  layer  thickness  keeps  growing.  As  a  result,  the 
normal  velocity  decreases.  The  maximum  normal  surface  velocity  occurs  at 
the  front  stagnation  point, while  the  minimum  is  located  near  the  separation 
point. 

Obviously,  droplet  vaporization  is  an  inherently  transient  process.  A 
summary  of  transient  droplet  life  history  is  presented  in  Figure  D-  4.  The 
results  indicate  that  both  transient  droplet  heating  and  reduction  in  Reynolds 
number  persist  during  most  of  the  droplet  lifetime  even  though  the  droplet 
has  reached  practically  wet-bulb  temperature  when  only  20%  of  the  mass  is 
vaporized.  These  two  transient  effects  constitute  the  major  sources  of  the 
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unsteady  behavior  of  this  problem. 
Overall  Transfer  Coefficients 


We  now  switch  our  attention  to  the  overall  characteristics  of  momentum, 
heat  and  mass  transfer  to  the  droplet.  Figure  D-  5  shows  the  variation  of 
drag  coefficients  of  droplets  initially  at  three  different  ambient  temperatures 
(Cases  1,  2,  and  3)  as  a  function  of  hydrodynamic  diffusion  time.  Also  shown 
in  this  figure  is  the  H-N-R  (Haywood,  Nafziger  and  Renksizbulut  [1989])  drag 
correlation  which  is  based  on  the  numerical  results  of  an  isolated  moving  drop 
vaporizing  in  its  own  fuel  vapor.  The  H-N-R  correlation  is  given  by, 
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CD(  1  +  BHJtlm)°-2  =  p— (1  +  0.2Re°m63);  10  <  Rem  <  300  (2) 


where  Bjjjiim  is  defined  to  account  for  droplet  heating.  The  gas  film  condi¬ 
tion  is  the  average  of  ambient  and  droplet  surface  conditions. 
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Due  to  the  sudden  injection  of  the  drop  into  a  uniform  flow-field,  the  time 
required  for  the  flow-field  to  relax  from  the  initially  impulsive  motion  can 
be  estimated  as  A tTeiax  ~  ~  0.4  diffusion  times.  During  the  initial  re¬ 

laxation  period  the  drag  coefficient  falls  rapidly.  Subsequently  Cd  increases 
as  a  result  of  a  reduction  in  the  Reynolds  number.  However,  for  the  higher 
ambient  temperature  case,  where  the  large  heat  flux  makes  droplet  vapor¬ 
ization  grow  fast  during  the  early  portion  of  its  lifetime,  the  drag  coefficient 
is  therefore  significantly  decreased  due  to  the  boundary  layer  blowing.  For 
all  cases,  Cd  tends  to  increase  with  time  during  the  final  portion  of  the  cal¬ 
culation,  which  implies  that  the  reduction  in  Reynolds  number  takes  control 
in  determining  drag  coefficients.  The  H-N-R  drag  correlation  does  not  seem 
applicable  in  our  cases.  The  deviation  can  be  as  much  as  20%  of  total  drag. 

Three  components  of  drag  coefficients  axe  shown  in  Figure  D-  6.  It  is 
clear  that  the  major  difference  in  total  drag  coefficients  for  the  three  different 
ambient  temperature  cases  comes  from  the  friction  drag  and  also  from  the 
pressure  drag  coefficients.  Higher  surface  blowing  can  reduce  the  friction 
drag  by  a  large  amount.  The  pressure  drag  coefficients  increase  steadily  as 
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a  result  of  the  reductions  in  upstream  velocity.  It  is  noteworthy  that  the 
thrust  drag  is  no  longer  negligible  for  the  high  transfer  number  case.  The 
contribution  from  the  —  VginVgigsin2d  term  may  account  for  more  than  a  10% 
in  reduction  of  the  total  drag.  This  term  physically  represents  the  thrust  force 
as  the  reaction  of  tangential  momentum  being  convected  radially  outwards 
and  serves  to  accelerate  the  droplet.  The  recoil  force  due  to  the  impact  of 
the  outward  flux,  Vg<n Vg>nsin20,  is  negligible. 

Figure  D-  7  shows  the  average  Nusselt  numbers  of  these  three  cases  and 
their  corresponding  H-N-R  Nusselt  number  correlations  which  can  be  ex¬ 
pressed  as 

Nuftlm(  1  +  B„Jllm)0-7  =  2  +  0.57  Rt]l2Pr)^m  (3) 

The  Nusselt  number  falls  during  the  initial  relaxation  period.  Further 
reduction  in  Nusselt  number  is  attributed  mainly  to  a  Reynolds  number  re¬ 
duction  and  an  increase  of  boundary  layer  thickness  due  to  surface  blowing. 
Again  our  numerical  results  do  not  agree  with  the  H-N-R  correlation.  How¬ 
ever,  the  discrepancy  becomes  smaller  during  the  final  part  of  the  calculation 
when  the  surface  temperature  approaches  the  wet-bulb  temperature. 

The  Sherwood  number,  which  represents  the  transient  dynamics  of  the 
mass  transfer,  is  shown  in  Figure  D-  8.  The  general  trend  of  variation  is 
very  similar  to  that  of  Nusselt  number.  The  relaxation  of  large  mass  fraction 
gradient  caused  by  the  sudden  injection  of  a  droplet  into  the  hot  gas-stream 
contributes  mainly  to  the  early  reduction  of  Sherwood  number.  A  more 
general  decrease  of  Sherwood  nun^r  is  then  observed,  which  may  be  due  to 
the  decay  of  Reynolds  number.  Our  numerical  value  is  lower  than  the  value 
predicted  by  the  correlation  of  H-N-R  which  is  given  by 

ShJilm{  1  +  BMJilm)07  =  2  +  0.87 Re^Sc1/^  (4) 

where  BM,jiim  =  ^{'lyj^ 

The  Parameter  Study 

Results  for  different  initial  droplet  temperatures  (300/f  vs  400A'  as  in 
Cases  1  and  8,  respectively)  are  shown  in  Figure  D-  9.  It  is  very  clear  that 
for  the  higher  initial  droplet  temperature  case,  the  surface  blowing  effect 
is  more  significant  than  that  of  the  lower  droplet  temperature  case  during  a 
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large  portion  of  the  droplet  lifetime.  Hence,  the  former  experiences  lower  drag 
than  the  latter.  Figure  D-  10  shows  the  portion  of  total  heat  flux  that  goes 
into  heating  of  the  drop  interior.  The  droplet  with  higher  initial  temperature 
would  spend  most  of  the  available  heat  energy  for  the  evaporation  process. 
On  the  contrary,  the  droplet  with  lower  initial  temperature  has  to  spend  most 
of  the  available  energy  to  heat  itself  first,  with  the  remaining  energy  for  va¬ 
porization.  The  heating  process  for  the  low  initial  droplet  temperature  case 
is  slower  and  will  persist  longer  than  for  the  high  initial  droplet  temperature 
case.  As  the  surface  temperatures  of  both  cases  approach  asymptotically  the 
same  wet-bulb  temperature  (i.e.  same  transfer  number),  the  drag  coefficient 
is  essentially  controlled  by  the  Reynolds  number.  Also  shown  in  Figure  D-  9 
is  the  solid  sphere  correlation  (see  Clift,  Grace,  and  Weber  (1978))  which  can 
be  grossly  inaccurate  due  to  the  high  mass  transfer  at  the  surface.  Since  the 
surface  blowing  effect  is  much  stronger  for  the  high  initial  droplet  tempera¬ 
ture  case,  the  Nusselt  and  Sherwood  numbers  are  smaller  than  those  for  the 
low  droplet  temperature  case. 

The  effect  of  variable  thermophysical  properties  is  very  important  for  the 
large  range  of  temperatures  considered  in  the  present  study.  A  calculation 
has  been  conducted  by  assuming  constant  properties  (except  for  gas-phase 
density)  that  are  based  upon  the  values  at  the  free  stream  (Cases  1  and  4). 
Since  the  global  behavior  is  evaluated  a  the  interface,  where  the  mixture 
gas-phase  composition  is  very  different  from  the  pure  gas  at  the  free  stream, 
the  constant  property  calculations  are  expected  to  predict  totally  different 
transfer  coefficients.  The  results  in  Figure  D-  9  show  that  the  constant 
property  case  can  overestimate  the  drag  coefficient  by  20%  or  more  compared 
to  the  variable  property  case.  The  discrepancy  is  attributed  to  the  change  in 
thermo-physical  properties  at  the  interface  and  the  change  in  the  flow  field 
caused  by  the  property  gradients. 

Figure  D-  9  shows  a  comparison  of  results  from  an  infinite  conductivity 
liquid-phase  model,  where  the  liquid-phase  temperature  is  assumed  to  be 
uniform  in  space  but  varying  with  time  (Cases  1  and  5).  The  no  slip  bound¬ 
ary  condition  is  also  applied  at  the  droplet  surface.  The  results  indicate  that 
the  infinite  conductivity  case  predicts  higher  drag  during  most  of  the  droplet 
lifetime.  For  the  infinite  conductivity  model,  the  liquid  phase  has  to  dis¬ 
tribute  all  the  available  heat  flux  to  the  droplet  interior.  Thus,  the  surface 
temperature  rises  slowly  resulting  in  a  very  weak  vaporization  rate,  as  well 
as  in  a  high  drag  coefficient.  The  transient  heating  period  is  relatively  short. 
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Because  of  the  larger  driving  temperature  potential,  the  droplet  eventually 
receives  more  heat  flux  than  that  calculated  from  our  model.  As  a  result, 
the  surface  temperature  does  eventually  increase  very  significantly,  and  the 
vaporization  grows  quickly  while  the  drag  coefficient  is  reduced  significantly. 
The  final  portion  of  drag  coefficient  is  primarily  determined  by  the  reduction 
of  Reynolds  number  as  already  illustrated  in  previous  cases. 

Calculations  for  three  different  initial  Reynolds  numbers  (Cases  1,  6,  and 
7)are  presented  in  Figure  D-  11.  It  is  observed  that  a  lower  initial  Reynolds 
number  results  in  a  higher  drag  coefficient  at  early  times  due  to  the  smaller 
convective  momentum  transport  which  usually  causes  the  lower  pressure  re¬ 
covery  at  the  rear  of  the  droplet.  It  is  also  noteworthy  that  the  three  drag 
coefficients  again  asymptotically  approach  a  certain  form  which  is  governed 
by  the  Reynolds  number  only.  This  is  because  the  same  effective  transfer 
number  has  been  reached  for  all  three  cases. 

The  final  parameter  study  is  concerned  with  the  fuel  volatility(Cases  1, 
9  and  10).  The  more  volatile  fuel  uses  a  large  fraction  of  the  energy  flux  to 
provide  the  latent  heat  for  vaporization  and  vaporizes  faster  than  the  heavier 
fuels.  The  less  volatile  fuel  spends  most  of  the  available  energy  to  heat  the 
droplet.  As  a  result,  the  surface  temperature  rises  very  slowly  and  the  effect 
of  surface  blowing  is  less  noticeable.  The  reduction  of  drag  due  to  surface 
blowing  is  more  pronounced  for  volatile  fuels  as  indicated  in  Figure  12. 

Correlations  for  the  Drag  Coefficient,  Nusselt  and  Sherwood  Numbers 

With  the  results  of  the  above  cases,  which  cover  a  wide  range  of  transfer 
numbers,  as  well  as  Reynolds  and  liquid  Peclet  numbers,  it  is  desirable  to 
obtain  the  correlations  of  drag  coefficient,  Nusselt  and  Sherwood  numbers  as 
a  function  of  some  important  parameters.  Since  these  transport  quantities 
are  determined  by  the  combined  influence  of  surface  blowing  and  vaporiza¬ 
tion,  the  internal  circulation,  and  the  unsteady  effects  related  to  the  droplet 
deceleration,  a  complete  correlation  which  covers  a  wide  range  of  all  param¬ 
eters  is  not  easy  to  achieve.  The  parameters  and  functional  form  selected  in 
this  work  are  the  same  as  in  the  H-N-R’s  correlation.  The  following  modi¬ 
fied  correlations  for  drag  coefficient,  Nusselt  and  Sherwood  numbers,  show  a 
good  agreement  (within  4%,  2%  and  6%  discrepancy,  respectively)  with  our 
numerical  results  as  shown  in  Figure  D-  5,  D-7  and  D-8. 
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Co(l  +  =  —(1  +  0.41ft£“») 
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Nu]ilm{\  +  BHJilm)0-7  =  2  +  0.46i?e^6  Prj{?m  (6) 

for  .6  <  finji/m  <  6.;  30  <  Rem  <  250 

The  correlation  of  Sherwood  number  is  approximated  by 

Shfilm(  1  +  BM,ftlm)0J  =  2  +  OAQRe^Sc1/^  (7) 

for  .3  <  BMjiim  <4.5 

The  major  difference  between  our  modified  drag  coefficient  correlation  and 
the  H-N-R  correlation  is  the  exponent  of  (1  +  Our  correlation  has  a 

higher  exponent  which  emphasizes  more  reduction  in  drag  coefficient  by  the 
effect  of  transfer  number. 


CALCULATION  OF  THE  TWO  DROPLETS  MOVING  IN  TANDEM 


In  the  present  research,  we  extend  the  study  of  Raju  and  Sirignano  [1988] 
to  include  the  effect  of  variable  properties.  The  single  droplet  code  is  modified 
with  the  addition  of  a  grid  generation  routine  and  a  treatment  to  the  read¬ 
justment  of  the  velocity  of  the  downstream  droplet,  as  well  as  the  variation 
of  spacing.  This  research  obviously  builds  upon  the  single  droplet  results. 

Table  D-3  and  Table  D-4  summarize  the  values  of  physical  parameters  in 
the  base  case  and  the  main  parameters  used  for  different  cases  examined  in 
the  present  research,  respectively.  Since  the  problem  involves  a  multidimen¬ 
sional  variables,  an  exhaustive  numerical  study  of  all  possible  combinations 
of  parameters  would  be  extremely  elaborate.  Therefore,  we  have  concen¬ 
trated  our  attention  on  a  detailed  study  of  the  interaction  effects  arising 
from  the  variation  of  Reynolds  number,  initial  diameter  ratio,  and  ratio  of 
initial  spacing  to  lead  droplet  diameter. 

The  results  are  presented  in  the  following  five  subsections.  The  first  three 
subsections  give  a  description  of  the  local,  as  well  as  overall,  behaviors  of  the 
interacting  droplets  for  the  base  case.  Most  of  the  results  show  that  the 
behavior  of  the  lead  droplet  is  qualitatively  in  agreement  with  that  of  iso¬ 
lated  droplet.  Hence  the  discussion  in  these  three  subsections  is  primarily 
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concentrated  on  the  behavior  of  the  downstream  droplet.  The  transient  vari¬ 
ation  of  the  spacings  for  different  initial  spacing,  initial  size  ratio,  and  initial 
Reynolds  numbers  is  given  in  the  fourth  subsection.  Subsection  5  presents  a 
correlation  for  the  drag  coefficients  of  the  lead  droplet  and  the  downstream 
droplet  from  the  constant  property  calculation. 

Results  for  the  Global  Flow  Field 

Figures  D-  13(a-g)  show  the  gas-phase  velocities,  as  well  as  the  contour 
plots  of  mass  fraction,  temperature,  vorticity  and  liquid-phase  streamlines, 
at  time  equals  2  when  the  spacing  is  reduced  from  8  to  6.  The  spacings  shown 
in  the  figures  are  nondimensionalized  with  respect  to  initial  radius  of  the  lead 
droplet.  The  gas-phase  velocities  plot  clearly  indicates  the  existence  of  recir¬ 
culation  zones  behind  the  droplets.  The  wake  length  for  the  leading  droplet 
is  about  one  droplet  diameter  which  is  in  close  agreement  with  the  single 
droplet  analysis.  There  is  only  a  small  recirculation  zone  formed  behind  the 
second  droplet.  Due  to  the  influence  of  the  recirculation  region  of  the  first 
droplet,  the  velocities  approaching  the  downstream  droplet  are  less  than  the 
velocities  upstream  of  the  lead  droplet.  As  a  result,  the  relative  velocity 
between  the  surrounding  gas  phase  and  the  downstream  droplet  is  consid¬ 
erably  smaller  than  that  of  the  lead  droplet.  It  can  be  expected  that  the 
downstream  droplet  will  receive  less  shear  stress  (less  momentum  transport 
from  gas  phase).  This  can  also  be  illustrated  from  the  diffusion  and  convec¬ 
tion  associated  with  the  gas-phase  vorticity  at  the  gas/liquid  interface.  The 
magnitude  of  vorticity  near  the  downstream  droplet  is  remarkably  smaller 
than  that  of  vorticity  at  the  region  close  to  the  lead  droplet,  as  shown  in  the 
vorticity  contours  of  Fig.  D-13b.  Also,  the  asymmetry  of  the  vorticity  dis¬ 
tribution  for  the  downstream  droplet  is  not  as  apparent  as  that  for  the  lead 
droplet.  This  means  that  the  downstream  droplet  receives  less  convective 
effect  from  the  gas  phase. 

In  Fig.  D-13c,  the  mass  fraction  contours  show  that  a  high  concentration 
gradient  occurs  near  the  lead  droplet  while  the  concentration  gradient  is 
lower  at  the  second  droplet.  The  fuel  vapor  generated  by  the  lead  droplet  is 
convected  downstream  and  alters  the  mixture  composition  of  the  surrounding 
environment  near  the  downstream  droplet.  The  relatively  rich  fuel  mixture 
is  usually  cold  enough  to  change  the  surrounding  thermal  environment  of 
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the  downstream  droplet  significantly.  The  thermal  boundary-layer  thickness 
of  the  downstream  droplet  is  larger  than  that  of  the  lead  droplet  as  shown 
in  the  isothermal  contour  plot  in  Fig.  D-13d.  Hence,  the  heat  and  mass 
transfer  rates  for  the  downstream  droplet  are  much  lower  than  those  of  the 
lead  droplet  due  to  the  presence  of  fuel  vapor. 

The  strength  of  the  stream  function  in  the  liquid  for  downstream  droplet 
is  almost  half  of  the  strength  for  the  lead  droplet.  This  may  be  primarily 
attributed  to  the  lower  shear  stress  acting  upon  the  downstream  droplet. 
The  liquid-phase  isotherms  of  Fig.  D-13g  indicate  that  a  thermal  boundary 
layer  exists  at  the  interface  for  the  lead  droplet  due  to  the  convection  effect. 
However,  for  the  downstream  droplet,  it  is  not  easy  to  identify  the  liquid 
thermal  boundary  layer.  The  heating  mechanism  for  the  downstream  droplet 
is  dominated  by  conduction.  The  smaller  magnitude  and  more  symmetric 
variation  are  two  main  characteristics  for  the  liquid-phase  vorticity  of  the 
downstream  droplet.  This  also  proves  that  the  external  flow  can  impose  only 
little  convective  influence  on  the  downstream  droplet. 

Local  Properties  Along  the  Droplet  Surface 

Figure  D-  14(a)  displays  surface  vorticity  at  different  times  as  a  function 
of  angular  position.  As  the  downstream  droplet  approaches  the  near  wake 
region  of  the  lead  droplet,  it  experiences  less  convection.  As  a  result,  the 
position  for  the  maximum  surface  vorticity  of  the  downstream  droplet  moves 
to  the  0  =  90°  plane.  Interestingly,  at  time  =  3,  the  surface  vorticity  at  the 
rear  portion  is  positive,  and  no  flow  separation  may  occur  there.  Since  the 
downstream  droplet  is  covered  by  the  wake  of  the  lead  droplet  at  time  =  3, 
the  downstream  droplet  is  in  a  very  low  Reynolds  number  flow  field  which 
makes  separation  impossible.  The  shear  stress  distribution  is  consistent  with 
these  conclusions  as  seen  in  Figure  D-  14(b).  The  negative  shear  stress  at  the 
stagnation  region  of  the  downstream  droplet  is  caused  by  the  recirculating 
flow  behind  the  lead  droplet. 

The  tangential  velocity,  as  seen  in  Figure  D-  14(c),  for  the  downstream 
droplet  is  much  smaller  than  that  of  the  lead  droplet  because  of  the  lower 
shear  stress  acting  on  the  surface. 

Figure  D-  14(d)  demonstrates  the  pressure  distribution.  Due  to  the  ac¬ 
tion  of  the  recirculating  wake,  the  pressure  at  the  stagnation  point  of  the 
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downstream  droplet  only  slightly  recovers  from  the  pressure  at  the  rear  stag¬ 
nation  point  of  the  lead  droplet.  As  the  spacing  decreases,  the  recovery  of 
the  pressure  at  the  stagnation  point  becomes  less.  It  is  expected  that  the 
pressure  drag  of  the  downstream  droplet  will  become  smaller  when  compared 
to  that  of  the  lead  droplet  as  the  droplets  approach  each  other. 

Figure  D-  14(e)  compares  the  local  Nusselt  number  variation  above  the 
droplet  surface  for  the  two  droplets.  The  decrease  of  the  Nusselt  number  at 
the  rear  part  of  the  lead  droplet  and  at  the  stagnation  region  of  the  down¬ 
stream  droplet  are  caused  by  the  approach  of  the  downstream  droplet.  As  the 
spacing  is  reduced,  the  increase  of  cold  mixture  density,  combined  with  the 
decrease  in  velocities  at  the  wake  region,  serve  to  retard  the  heat  exchange. 
The  increase  of  Nusselt  number  at  the  rear  part  of  the  downstream  droplet 
is  caused  by  the  hot  ambient  stream  entrained  by  the  recirculating  flow.  A 
similar  trend  can  be  observed  for  the  local  Sherwood  number  as  shown  in 
Figure  D-  14(f).  Figure  D-  14(g)  gives  the  surface  temperature  distribution. 
The  highly  nonuniform  distribution  of  the  surface  temperature  for  the  down¬ 
stream  droplet  is  due  to  the  slow  circulation  which  cannot  distribute  energy 
efficiently  in  the  streamwise  direction.  As  a  result,  the  distribution  of  surface 
temperature  is  similar  to  that  of  the  local  Nusselt  number.  The  heat  flux  is 
diffused  in  the  radial  direction  within  the  downstream  droplet.  It  may  cause 
a  small  cold  region  at  the  front  and  a  large  hot  region  at  the  rear  of  the 
droplet. 

Overall  Behaviors  of  the  Two  Droplets 

In  this  subsection,  the  drag  coefficient,  Nusselt  and  Sherwood  numbers, 
and  average  volumetric  temperature  for  the  lead  droplet  and  the  downstream 
droplet  are  presented.  The  results  for  an  isolated  droplet  in  the  same  convec¬ 
tive  flow  field  are  also  presented  for  comparison.  Figure  D-  15  shows  that  the 
drag  coefficient  of  the  lead  droplet  drops  about  6  %  from  its  isolated- droplet 
value  due  to  the  interaction  with  the  downstream  droplet.  The  discrepancy 
increases  as  the  droplet  spacing  decreases.  The  downstream  droplet  expe¬ 
riences  much  lower  drag  because  of  the  wake  effect  from  the  lead  droplet. 
The  three  drag  components  are  shown  in  Figure  D-  16.  The  difference  in 
friction  drag  seems  to  contribute  most  of  the  drag  difference  between  the 
lead  droplet  and  the  isolated  droplet.  The  sharp  drop  of  the  total  drag  at 
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the  time  when  the  two  drops  are  about  to  collide  is  mainly  attributed  to 
the  reduction  in  pressure  drag.  Figure  D-  17  and  Figure  D-  18  present  the 
overall  Nusselt  numbers,  Nu,  and  average  volumetric  temperature,  Tv,  re¬ 
spectively.  The  difference  in  both  Nu  and  Tv  between  two  droplets  indicate 
that  wake  effects  tend  to  reduce  significantly  the  heat  transport  to  the  down¬ 
stream  drop.  For  the  lead  droplet,  the  influence  of  the  downstream  droplet  is 
not  negligible  even  when  the  spacing  is  about  four  droplet  diameters.  Other 
results  for  the  cases  of  separating  drops  show  that  the  Nusselt  numbers  for 
the  lead  and  downstream  droplets  behave  as  for  the  isolated  droplet  during 
the  final  period  of  calculation.  The  Sherwood  number  of  the  downstream 
droplet,  shown  in  Figure  D-  19,  indicates  a  longer  relaxation  period  for  the 
mass  exchange  process,  which  is  caused  by  the  upstream  droplet,  throughout 
our  calculation. 


The  Transient  Variation  of  Spacing 


The  determination  of  whether  two  droplets  are  going  to  collide  or  separate 
after  they  are  suddenly  introduced  into  a  combustor  is  of  great  interest.  It 
can  be  easily  shown  that  the  droplet  deceleration  is  proportional  to  its  drag 
coefficient  and  inversely  proportional  to  its  radius.  The  change  of  spacing  is 
proportional  to  the  relative  deceleration  between  two  droplets.  This  can  be 
expressed  as  follows: 

AS  (differential  change  of  spacing)  oc  (Ad  2  —  Ad\)*  (A Tng)2 

where  Ad  is  the  droplet  deceleration  and  A r//5  is  a  computational  time  step 
in  gas  phase  hydrodynamic  diffusion  time  scale.  The  subscripts  1,2  corre¬ 
spond  to  the  lead  and  downstream  droplet,  respectively.  Positive  AS  means 
that  droplet  separation  occurs.  Negative  AS  means  that  droplet  coalescence 
prevails.  In  general,  Cd,i  is  larger  than  Cd, 2-  If  R2  is  larger  than  or  equal 
to  R\ ,  droplets  will  collide.  On  the  contrary,  if  R2  is  much  smaller  than  Rlt 
droplet  separation  is  most  likely  to  occur. 

The  transient  variation  of  center-to-center-spacings  for  the  cases  of  dif¬ 
ferent  initial  spacing  is  depicted  in  Figure  D-”20.  Droplet  collision  appears 
to  be  likely  for  initially  equal-sized  droplets,  as  long  as  the  initial  spacing  is 
sufficiently  small  for  the  wake  effects  to  be  important.  The  drag  coefficients 
of  two  droplets  for  initial  spacing  of  4,  8  and  16  are  shown  in  Figure  D-  21. 
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The  smaller  the  initial  spacing,  the  more  the  drag  coefficient  deviates  from  its 
isolated  value  for  the  lead  droplet  and  the  more  the  drag  coefficient  decreases 
for  the  downstream  droplet.  Figure  D-  22  demonstrates  how  two  droplets 
will  separate  or  collide  for  the  cases  of  different  initial  size  ratio.  The  first 
observation  is  that  the  larger  the  difference  in  droplet  size,  the  faster  they  ap¬ 
proach  each  other(or  separate).  For  the  case  of  R2  =  0.75,  droplet  separation 
occurs  at  a  very  early  time,  but  eventually  droplet  coalescence  prevails.  The 
results  of  Raju  and  Sirignano  [1988]  indicate  that  there  exists  a  bifurcation 
point  depending  upon  the  critical  ratio  of  R2io/R'ifi  below  which  droplet  co¬ 
alescence  becomes  unlikely.  They  predict  that  R2  =  0.75  is  very  close  to  the 
bifurcation  point,  so  that  the  droplets  spend  more  time  under  the  influence 
of  each  other.  Figure  D-  23  displays  the  drag  coefficients  of  two  droplets 
for  the  cases  of  R2  <  1.  The  drag  coefficients  of  the  lead  droplets  for  the 
cases  of  smaller  downstream  droplet  are  all  collapsed  together.  This  means 
at  this  spacing  the  lead  droplet  receives  little  influence  from  the  downstream 
droplet.  The  drag  coefficients  of  the  downstream  droplets  for  the  cases  of 
droplet  separation  tend  to  increase  with  the  reduction  in  Reynolds  number. 
This  is  also  the  typical  characteristic  of  an  isolated  droplet. 

The  influence  of  the  initial  Reynolds  number  for  the  case  of  droplet  co¬ 
alescence  is  shown  in  Figure  D-  24.  An  increase  in  initial  Reynolds  number 
seems  to  increase  the  approach  speed.  Figure  D-  25  indicates  that  by  de¬ 
creasing  the  initial  Reynolds  number  for  the  case  with  droplet  size  ratio  at 
the  critical  point,  the  droplet-motion  pattern  can  be  changed  from  coales¬ 
cence  to  separation.  However,  for  the  case  of  small  initial  spacing,  the  initial 
Reynolds  number  does  not  play  a  role  as  shown  by  the  results  in  Figure  D- 
26.  From  the  results  shown  above,  the  initial  spacing  and  initial  Reynolds 
number  are  not  as  important  as  the  droplet  size  ratio  in  the  determination 
of  droplet  separation  or  coalescence. 

Since  the  droplet  interaction  effects  are  determined  by  the  combined  in¬ 
fluence  of  surface  blowing  and  vaporization,  the  internal  circulation,  and  the 
unsteady  effects  related  to  the  droplet  deceleration,  a  parameter  study  must 
cover  a  wide  range  of  initial  Reynolds  number,  size  ratio,  initial  spacing, 
transfer  number,  etc.,  in  order  to  obtain  a  meaningful  set  of  correlations  for 
the  transfer  coefficients.  The  correlations  of  the  drag  coefficient  and  Nus- 
selt  and  Sherwood  numbers  for  the  lead  and  downstream  droplet  are  still 
under  investigation.  Raju  and  Sirignano  [1988]  provided  the  following  cor¬ 
relation  for  the  drag  coefficients  of  two  droplets  from  the  constant  property 
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calculation. 


cD3  =  (^1)(^)CD.,  - 1.22 

'oo,2  ^1  '00,2 


'3,0 


(8) 


where  V^,,-  is  the  relative  velocity  between  far  upstream  gas  phase  and  the 
2-th  droplet. 

The  numerical  results  satisfy  the  above  equation  within  12  percent  ac¬ 
curacy  for  the  range  of  initial  spacings  between  4-30  and  initial  Reynolds 
numbers  between  50-200. 
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Table  D-l  VALUES  OF  PHYSICAL  PARAMETERS  USED  IN  THE  ISOLATED  DROPLET 

BASE  CASE  COMPUTATION 


Parameter 

Value 

Initial  Reynolds  number,  gas  phase  Re  =  2a'0U'!X>0p'g <<x>/ p'^ 

100.0 

Relative  velocity  of  drop,  m/s 

25.0 

Free  stream  temperature,  [K] 

1250.0 

Combustor  pressure,  atm 

10.0 

Prandtl  number,  gas  phase 

0.74 

Prandtl  number,  liquid  phase 

8.59 

Schmidt  number,  gas  phase 

2.36 

Molecular  weight,  oxidizer,  Kg/Kmol 

29.0 

Molecular  weight,  fuel,  n-octane,  Kg/Kmol 

114.2 

Droplet  initial  temperature,  [K] 

300.0 

Viscosity  ratio,  p'J  p'g<oc 

10.49 

Density  ratio,  p'i/p'g<00 

251.93 

Specific  heat  at  constant  pressure  ratio, Cp|0/Cp^  ^ 

1.87 

Latent  heat  /  Specific  heat  of  liquid,  [I<] 

135.95 

Table  D-2  MAIN  PARAMETERS  CONSIDERED  IN  EACH  CASE  OF  THE  ISOLATED 

DROPLET  STUDY 


|g 

Ambient 
Temperature  [I<] 

Initial  Droplet 
Temperature  [K] 

Initial 
Reynolds  # 

Fuel 

Type 

Droplet 

Heating  Model 

Physical 

Property 

i. 

1 

100 

n-octane 

two-dimensional 

variable 

2. 

100 

n-octane 

two-dimensional 

variable 

3. 

300 

100 

n-octane 

two-dimensional 

variable 

4. 

1250 

300 

100 

n-octane 

two-dimensional 

constant 

5. 

1250 

300 

100 

n-octane 

oo  conductivity 

variable 

6. 

1250 

300 

150 

n-octane 

two-dimensional 

variable 

7. 

1250 

300 

50 

n-octane 

two-dimensional 

variable 

8. 

1250 

400 

100 

n-octane 

two-dimensional 

variable 

9. 

1250 

300 

100 

n-decane 

two-dimensional 

variable 

10. 

1250 

300 

100 

n-hexane 

two-dimensional 

variable 

51 


Table  D-3  VALUES  OF  PHYSICAL  PARAMETERS  USED  IN  THE  TWO  DROPLETS 

BASELINE  COMPUTATION 


Parameter 

Value  || 

Initial  Reynolds  number,  gas  phase  Re  =  2a'0U'oo  0p'g tZJ 

100.0  II 

Relative  velocity  of  drop,  m/s 

2-5. G 

Free  stream  temperature,  [K] 

1000.0 

Combustor  pressure,  atm 

10.0 

Prandtl  number,  gas  phase 

0  714 

Prandtl  number,  liquid  phase 

14.92 

Schmidt  number,  gas  phase 

3.26 

Molecular  weight,  oxidizer,  Kg/Kmol 

29.0 

Molecular  weight,  fuel,  n-octane,  Kg/Kmol 

142.28 

Droplet  initial  temperature,  [K] 

300.0 

Viscosity  ratio,  p'J p'goo 

21.44 

Density  ratio,  p\!  p'g^ 

209.19 

Specific  heat  at  constant  pressure  ratio, Cp\  0/Cp'g  ^ 

1.94 

Latent  heat  /  Specific  heat  of  liquid,  [I<] 

126.28 

Table  D-4  MAIN  PARAMETERS  CONSIDERED  IN  EACH  CASE  OF  THE  TWO  DROPLETS 

COMPUTATION 


Case 

# 

Ambient 
Temperature  [K] 

Initial  Droplet 
Temperature  [K] 

Fuel 

Type 

Initial 
Reynolds  # 

d'/(r'0) 

spacing 

Physical 

Property 

1. 

1000 

300 

n-decane 

100 

1 

8 

variable 

2. 

1000 

300 

n-decane 

100 

1 

16 

variable 

3. 

1000 

300 

n-decane 

100 

1 

32 

variable 

4. 

1000 

300 

n-decane 

100 

1 

4 

variable 

5. 

1000 

300 

n-decane 

50 

1 

8 

variable 

6. 

1000 

300 

n-decane 

25 

1 

8 

variable 

7. 

1000 

300 

n-decane 

10 

1 

8 

variable 

8. 

1000 

300 

n-decane 

125 

1 

8 

variable 

9. 

1000 

300 

n-decane 

100 

H 

8 

variable 

10. 

1000 

300 

n-decane 

100 

EB 

8 

variable 

11. 

1000 

300 

n-decane 

100 

0.4 

8 

variable 

12. 

1000 

300 

n-decane 

100 

0.2 

8 

variable 

13. 

1000 

300 

n-decane 

100 

1.2 

8 

variable 

14. 

1000 

300 

n-decane 

100 

0.7 

4 

variable 

15. 

1000 

300 

n-decane 

100 

0.5 

4 

variable 

16. 

1000 

300 

n-decane 

50 

0.7 

4 

variable 

17. 

1000 

300 

n-decane 

50 

0.5 

4 

variable 

18. 

1000 

300 

n-decane 

100 

0.65 

8 

variable 

19. 

1000 

300 

n-decane 

30 

OB 

8 

variable 
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SHEAR  shear 

AIR  I  LIQUID  AIR  II 


Figure  A-l:  Two-Dimensional  Parallel  Flow  Plane  Sheet  Atomizer 
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Figure  A-4:  “Stretched  Streamwise  Ligament  Breakup”:  Water  Sheet  at 
1  m/s  with  20 m/a  Co-Flow  Air:  (a)  Face  View,  (b)  Side  View. 
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Figure  A-5:  Drop  Size  Distribution:  (a)  Cellular  Breakup,  (b)  Stretched 
Streamwise  Ligament  Breakup 
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Figure  A-6:  Effect  of  Relative  Velocity  on:  (a)  Sheet  Breakup  Length,  (b) 
Sheet  Breakup  Time. 


Figure  A-7:  Effect  of  Liquid  Properties  on  “Cellular  Breakup”:  (a)  Face 
View,  Water  and  Ethanol  (Case  13)  (b)  Face  View,  Water  and  Ethanol  (Case 
26),  (c)  Face  View,  Water  and  Ethylene  Glycol  (Case  14). 
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Figure  A-8:  Effect  of  Liquid  Properties  on  ” Stretched  Streamwise  Ligament 
Breakup”  (Case  22):  (a)  Water,  (b)  Ethanol,  (c)  Ethylene  Glycol. 
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Figure  B-2:  Evolution  of  the  interface  shape  (a)  and  the  interface  vorticitv 
(b)  for  p  =  1 .  IF  =  0.1 
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Figure  B-4:  Maximum  vortex  strength  vs.  time  for  p  =  1  and  selected  values 
of  W. 


Figure  B-5:  Maximum  vortex  strength  vs.  time  for  W  =  2/3  and  selected 
values  of  p. 


Figure  B-6:  Shear  instability  of  a  planar  fluid  sheet. 
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Figure  B-7:  Sheet,  evolut  ion  for  a  thickness-to-wavelcngtli  ratio  of  2.  Density 
ratio  =  1.  Surface  tension  =  0.  (shiftcd-symmetric  mode). 
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Figure  B-9:  Sheet  evolution  for  a  thickness- to- wavelength  ratio  of  1/4.  Den¬ 
sity  ratio  =  1.  No  surface  tension,  (symmetric  mode) 
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Figure  C-l:  Evolution  of  the  interface  for  Re  =  50 
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Figure  C-4:  Species  Concentration  PDF  constructed  from  the  profiles  in 
Figure  C-2. 
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Figure  C-5:  PDF  for  Re  =  10  and  Sc  =  1. 
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Figure  C-6:  PDF  for  Re  —  100  and  Sc  =  1. 


Figure  C-9:  Plot  of  the  relative  difference  R  between  the  PDF  on  the  low- 
speed  and  high-speed  streams,  a)  Re  =  505c  =  1000  b)  Re  =  105c  =  10  c) 
Re  =  1005c  =  1 


!  C-12:  Mixedness 
the  mixing  plane. 
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Figure  C-13:  Comparison  of  the  overalll  mixedness  parameter  between 
mixing  layer  and  an  isolated  vortex. 


contour  interval:  2.54E-02  from: -5.23E-02  to:  4.56E-01  t  —  1.25 


contour  interval:  2.97E-02  from: -5.00E-02  to:  5.44E-01  * 


contour  interval:  3.24E-02  from: -5.00E-02  to:  5.97E-01  <  —  1-75 


Figure  C-14:  Streamfunction  contours  for  the  simulated  Mixing  Layer. 
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Figure  C-15:  Rollup  of  an  initially  sinusoidal  disturbance,  (a)  Rediscretiza¬ 
tion  and  addition  of  vortices,  (b)  Combination  and  addition  of  Vortices. 
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Figure  C-16:  Rollup  of  pairing  of  two  vortical  structures  induced  by  the 
presence  of  the  first  subh&rmonic. 
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Figure  C-17:  Rollup  and  grouping  of  three  vortical  structures  induced  by  the 
presence  of  the  second  subharmonic. 
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Figure  C-18:  Evolution  of  the  vorticity  distribution  along  the  interface  during 
the  grouping  of  the  three  vortical  structures. 
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Figure  C-19:  Amplitude  of  the  disturbance  as  a  function  of  time  for  the  fun¬ 
damental  disturbance  alone  (0),  the  fundamental  plus  the  first  subharmonic 
(1),  and  the  fundamental  plus  the  second  subharmonic  (2). 
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Figure  C-21:  Evolution  of  Mass  Fraction  Contours  with  four  structures. 


GAS  PHASE  VELOCITY  VECTORS 


Time  =  25,00 
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Ambient  Temperature =  2250  ^ 
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Initial  Reynolds  Number  =  200 

Figure  D-l(a)  Instantaneous  gas  phase  velocity  at  time  =  25. 
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Figure  D-l(b)  Mass  fraction  contour  plot  at  time  =  25. 
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Figure  D-l(c)  Temperature  contour  plot  at  time  =  25. 
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Figure  D-l(d)  Stream  function  of  liquid  phase  at  time  =  25. 
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Figure  D-l(e)  Vorticity  contour  plot  at  time  =  25. 
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Figure  D-2(&)  Transient  history  of  droplet  heating.  Time  =  0.5 


LIQUID  PHASE  ISOTHERMS 

Contour  Interval:  5.76E+O0  Min:  3.82E+02  Max:  4.68E+02 


Reynolds  Number-  59,90 


Ambient  Temperature =  1250  K 

Initial  Droplet  Temperature =  300  K 
Initial  Reynolds  Number  =  100 
Fuel :  n-octane 


Figure  D-2(b)  Transient  history  of  droplet  heating.  Time  =  10.0 
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Figure  D-2(c)  Transient  history  of  droplet  heating.  Time  =  25.0 
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Figure  D*3(a)  Surface  shear  stress  distribution  at  different  times. 
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Figure  D-3(b)  Surface  tangential  velocity  distribution  at  different  times. 
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Initial  Reynolds  Number  s  100. 


Figure  D-3(c)  Surface  nondimensional  pressure  distribution  at  different  times. 
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_ t  =  1.0,  Reft)  =  93.6,  S.P.  at  156.3,  Local  min.  at  159.5 

_ t  =  3.0,  Reft)  =  83.7,  S.P.  at  180.0,  Local  min.  at  180.0 

_ t  =  6.0,  Reft)  =  72.4,  S.P.  at  180.0,  Local  min.  at  180.0 

_  t=  9.0,  Reft)  =  623,  S.F. at  180.0,  Local  min. at  180.0 

_ t  =  12.0,  Reft)  =  54.2,  S.P.  at  167.9,  Local  min.  at  165.6 

_ t=  15.9,  Reft)  =  43.9,  S.P.  at  157.0,  Local  min.  at  153.2 

_ t=  20.0,  Reft)  =  34.3,  S.P.  at  149.2,  Local  min.  at  147.3 

_ t  =  25.0,  Reft)  =  24.1,  S.P.  at  141.3,  Local  min.  at  146.9 
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Figure  D-3(d)  Local  Nusselt  number  distribution  at  different  time*. 
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ANGULAR  POSITION  (DEGREE) 

Fual :  n-octana 

Ambiant  Tampa  ratura  81250.  K 
Initial  Droplat  Tamparatura  m  300.  K 
Initial  Reynolds  Number  m  100. 

Figure  D*3(f)  Surface  normal  velocity  distribution  at  different  times. 


DRAG  COEFFICIENTS 


_ 1.  Tg =1250  K,  11=300  K,  n-octane  fuel,  Numerical  Solution. 

_ 2.  To =1000  K,  TfcXO  K,  n-octane  fuel,  Numerical  Solution. 

_ 3.  Tfl=  800  K,  Tfe300  K,  n-octane  Juel,  Numerical  Solution. 

_ 4.  Tq =1250  K,  Renkslzbulut-Haywood  Correlation 

_ 5.  To =1000  K,  Renkslzbulut-Haywood  Correlation 

_ 6.  Tq=  800  K,  Renkslzbulut-Haywood  Correlation 

_ 7.  To =1250  K,  Modified  Correlation 

_ 8.  To =1000  K,  Modified  Correlation 

_ 9.  Tfl=  800  K,  Modified  Correlation 


GAS  HYDRODYNAMIC  DIFFUSION  TIME  SCALE 


Figure  D-5  Time  variation  of  drag  coefficients  for  different  ambient  temperatures  as  well  as  com¬ 
parison  with  H-N-R  correlation  and  modified  correlation. 


DRAG  COEFFICIENT  COMPONENTS 


_ 1.  Tg=l250  K,  Tl=300  K,  n-octane  fust,  Pressure  Drag. 

_ 2.  Tfl=1000  K,  T1=30GL  K,  n-octane  fuel,  Pressure  Drag. 

_ 3.  Tg=  800  K,  Tl=300  K,  n-octane  .fuel,  Prassura  Drag. 

_ 4.  Tg=12S0  K,  Tl=300  K,  n-octane  fuel,  Friction  Drag. 

_ 5.  Tgsiooo  K,  Ti=300  K,  n-octane  fuel,  Friction  Drag. 

_ 6.  Tg=  800  K,  Ti=300  K,  n-octana  fual,  Friction  Drag. 

_ 7.  Tg=1250  K,  Ti=300  K,  n-octane  fual,  Thrust  Drag. 

_ 8.  Tg=lOOO  K,  T)s300  K,  n-octana  fual,  Thrust  Drag. 

_ 8.  Tg=  800  K,  Tis300  K,  n-octana  fual.  Thrust  Drag. 


GAS  HYDRODYNAMIC  DIFFUSION  TIME  SCALE 


Figure  D-6  Time  variation  of  three  drag  components  for  different  ambient  temperatures. 


NUSSELT  NUMBER 


_  i.  Tg=1250  K,  Tl=300  K,  n-octane  fuel,  Numerical  Solution. 
.  2.  TgslOOO  K,  Tl=300  JC,  n-octana  fuel,  Numerical  Solution. 
_  3.  Tg=  800  K,  Tis300  K,  n-octane  fuel,  Numerical  Solution. 
_  4.  Tg=1250  K,  Renksizbulut-Haywood  Correlation 
_  5.  TgslOOO  K,  Renksizbulut-Haywood  Correlation 
.  6.  Tg=  800  K,  Renksizbulut-Haywood  Correlation 
.  7.  Tg=12S0  K,  Modified  Correlation 
.  8.  Tg=1000  K,  Modified  Correlation 
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GAS  HYDRODYNAMIC  DIFFUSION  TIME  SCALE 


Figure  D-7  Time  variation  of  average  Nusselt  numbers  for  different  ambient  temperatures  as  well 
as  comparison  with  H-N-R  correlation  and  modified  correlation. 


SHERWOOD  NUMBER 


_ 1.  Tg=l250  K,  Tl=300  K,  n-octane  fuel,  Numerical  Solution. 

_ 2.  Tgsiooo  K,  TIsSOO  K,  n-octane  fual,  Numerical  Solution. 

_ _ 3.  Tg=  800  K,  11=300  K,  n-octane  .fuel,  Numerical  Solution. 

_ 4.  Tg=1250  K,  Renkslzbulut-Haywood  Correlation 

_ 5.  TgslOOO  K,  Renkslzbulut-Haywood  Correlation 

_ 6.  Tg=  800  K,  Renkslzbulut-Haywood  Correlation 

_ 7.  Tg=1250  K,  Modified  Correlation 

_ 8.  Tg=1000  K,  Modified  Correlation 

_ 9.  Tg=  800  K,  Modified  Correlation 


GAS  HYDRODYNAMIC  DIFFUSION  TIME  SCALE 


Figure  D-8  Time  variation  of  average  Sherwood  Dumber*  for  different  ambient  temperatures  as  well 
as  comparison  with  H-N-R  correlation  and  modified  correlation. 


DRAG  COEFFICIENTS 


_  1.  Tis300  K,  Variable  Property,  Axlsymmetric  Modal. 

.  2.  Tis400  K,  Variable  Property,  Axlsymmetric  Model. 

_  3.  T1=300  K,  Constant  Property,  Axlsymmetric  Model. 

_  4.  Tls300  K,  Variable  Property,  Infinite  Conductivity  Model. 
.  5.  Solid  Sphere  Drag  Correlation. 


INSTANTANEOUS  REYNOLDS  NUMBER 


Figure  D-9  Time  variation  of  drag  coefficients  for  the  cases  of  different  initial  droplet  temperatures, 
different  droplet  beating  model  and  constant  property  assumption. 


ATING 


1.  Tg=1250  K,  T1=300  K,  n-octane  fuel,  Numerical  Solution. 

2.  Tg=1250  K,  TI=400  K,  n-octane  fuel,  Numerical  Solution. 
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GAS  HYDRODYNAMIC  DIFFUSION  TIME  SCALE 


Figure  D-10  Time  variation  of  the  rate  of  droplet  heating  for  different  initial  droplet  temperatures 


DRAG  COEFFICIENTS 


_ 1.  Tg=1250  K,  Tl=300  K,  Re=150,  n-octane  fuel,  Numerical  Solution. 

_ 2.  Tg=1250  K,  Tl=300  K,  RoslOO,  n-oetane  fuol,  Numerical  Solution. 

_ 3.  To =1250  K,  Tl=300  K,  Re=50,  n-octane  fuel,  Numerical  Solution. 
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INSTANTANEOUS  REYNOLDS  NUMBER 


Figure  D-ll  Time  narration  of  drag  coefficient.  for  different  initial  Reynold,  number.. 


DRAG  COEFFICIENTS 


.  1.  Tfl=1250  K,  11=300  K,  n-hexane  fuel,  Numerical  Solution. 
..  2.  Tg=1250  K,  T1=300  K,  n-octane  fuel,  Numerical  Solution. 
.  3.  Tg=1250  K,  11=300  K,  n-decane  fuel,  Numerical  Solution. 
.  4.  Solid  Sphere  Drag  Correlation. 


INSTANTANEOUS  REYNOLDS  NUMBER 


Figure  D-12  Time  variation  of  drag  coefficients  for  different  fuels. 


GAS  PHASE  VELOCITY  VECTORS 


At  time  =  2.00 

Re  #(1)  =  85.95  ,  Re  #(2)  =  P0.35 
7?1  =  7.00  ,  R2  =  7.00  ,  Spacing  =  5.00 
Ambient  Temperature-  1000  K 
Initial  Droplet  Temperature-  300 
Initial  Reynolds  Number  -  100 
Fuel :  n-decane 
Initial  Rl=l,  R2-1 ,  Spacing=8 


Figure  D-13(a)  Instantaneous  gas  phase  velocity  at  time  =  2. 


GAS  AND  LIQUID  PHASE  VORTICITY 

Contour  Interval:  3.10E-01  Min:  -5.37E+00  Max:  2.03E-01 


Time  =  2.00 

Re  #(1)  =  85.95  ,  Re  #(2)  =  90.35 
R1  =  1.00  ,  R2  -  1.00  ,  Spacing  = 

Ambient  Temperature =  7000  F 
Initial  Droplet  Temperature-  300  K 
Initial  Reynolds  Number  =  100 
File/ :  n-decane 
Initial  Rl-1 ,  F2=7,  Spacing=8 


Figure  D- 13(b)  Vortidty  contour  plot  at  time  =  2. 


MASS  FRACTION 


Time  =  2.00 

Re  #(1)  =  85.95  ,  Re  %(2)  =  90.35 
1W  =  1.00  ,  R2  =  1.00  ,  Spacing  =  5.00 
Ambient  Temperature-  1000  K 
Initial  Droplet  Temperature =  300  ]£ 

Initial  Reynolds  Number  -  100 
Fuel :  n-decane 
Initial  Rl=l,  R2=l,  Spacing=8 


Figure  D-13(c)  Mass  fraction  contour  plot  at  time  =  2. 


Time  —  2.00 

Re  mi)  =  85.95  ,  Re  *(2)  =  90.35 
R1  =  7.00  ,  R2  =  7.00  ,  Spacing  =  0.00 
Ambient  Temperature =  7000  AT 

Initial  Droplet  Temperature -  300  K 
Initial  Reynolds  Number  =  700 
Fwe/ ;  n-decane 
Initial  Rl=l,  R2-1,  Spacing=8 


Figure  D-13(d)  Temperature  contour  plot  at  time  =  2. 


laid 


ST  DROPLET  LIQUID  PHASE  STREAM  FUNCTION 


Contour  Interval:  1.68E-01  Min:  O.OOE+OOMax:  3.36E+00 


1.68E-01. 


SECOND  DROPLET  LIQUID  PHASE  STREAM  FUNCTION 


Contour  Interval:  8.25E-02  Min:  O.OOE+OOMax:  1.65E+00 


9.90E-01 


8.25E-02. 


Re#(l)=  85.95  ,  Re  %(2)  =  90.35 
R1  =  1.00  ,  R2  =  1.00  ,  Spacing  =  6.00 
Ambient  Temperature^  1000  K 
Initial  Droplet  Temperature^  300 
Initial  Reynolds  Number  =  100 
Fuel :  n-decane 
Initial  Rl=l,  R2=l,  Svacine=8 

Figure  D- 13(e)  Stream  function  of  liquid  phase  at  time  =  2. 


FIRST  DROPLET  LIQUID  PHASE  VORTICITY 

Contour  Interval:  4.20E-02  Min: -7.42E-01  Max:  1.43E-02 


SECOND  DROPLET  LIQUID  PHASE  VORTICITY 

Contour  Interval:  1.66E-02  Min: -2.80E-01  Max:  1.88E-02 
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17.5  18.5 

Time  =  2.00 

Re  #(1)  =  85.95  ,  Re  #(2)  m  9035 
R1  =  1.00  ,  R2  =  1.00  ,  Spacing  =  0.00 
Ambient  Temperatures  1000  K 
Initial  Droplet  Temperature =  300 
Initial  Reynolds  Number^  100 
Fuel :  n-decane 
Initial  Rlsl,  R2=l,  Spacing s8 

Figure  D- 13(f)  Detailed  vorticity  contour  of  liquid  phase  at  time  =  2. 


FIRST  DROPLET  LIQUID  PHASE  ISOTHERMS 


Contour  Interval:  4.56E+00  Min:  3.04E+02Max:  3.77E+02 

l  n  _ _ _ — 


SECOND  DROPLET  LIQUID  PHASE  ISOTHERMS 


Contour  Interval:  7.64E+00  Min:  3.03E+02  Max:  4.25E+02 

l.Q_ _ — - 


Time  =  2.00 

Re  #(1)  =  85.95  ,  Re  #(2)  =  90.35 


R1  =  1.00  ,  R2  =  1.00  ,  Spacing  =  6.00 
Ambient  Temperature^  1000  K 
Initial  Droplet  Temperature-  300  K 
Initial  Reynolds  Number  -  100 
Fuel :  n  decane 
Initial  Rl-1 ,  R2=l,  Spacing=8 

Figure  D- 13(g)  Detailed  temperature  contour  of  liquid  phase  at  time  =  2. 


ANGULAR  POSITION  (DEGREE) 

Fuel :  N-DECANE 
Ambient  Temperature  *1000.  K 
Initial  Droplet  Temperature  *  300.  K 
Initial  Reynolds  Number  ■  100. 

Initial R1  * IjOO  ,  R2-1.00  ,  D«8.0 

Figure  D-14(a)  Surface  vortidty  distribution  at  different  times  for  the  lead  droplet  and  toe  oown 
stream  droplet. 


_  Lead  droplet,  t=  1.0,  Re(t)  =  92.2,  S.P.  at  141.4,  Stress=0  at  122.1 

_ Lead  droplet,  t*  2.0,  Re(t)  =  85.9,  S.P.  at  148.9,  Stress=0  at  122.3 

_ Lead  droplet,  t»  3.0,  Re(t)  =  80.7,  S.P.  at  161.2,  Stress=0  at  118.3 

_ Downstream  droplet,  t  =  1.0,  Re(t)  =  94.2,  S.P.  at  144.7,  StresssO  at  133.9 

_ Downstream  droplet,  t  =  2.0,  Re(t)  =  90.3,  S.P.  at  157.6,  StresssO  at  149.1 

_ Downstream  droplet,  t  s  3.0,  Re(t)s87.1,S.P.at  12.0,  Stress=0  at  27.7 


ANGULAR  POSITION  (DEGREE) 

Fuel :  N-DECANE 
Ambient  Temperature  *1000.  K 
Initial  Droplet  Temperature  *  300.  K 
Initial  Reynolds  Number  ■  100. 

Initial R1  ■  1j00  ,  R2.1JJ0  ,  D-8.0 

Figure  D.  14(b)  Surface  .beer  rtrere  dirtributiou  .t  diterent  time,  fcr  the  W  droplet  tod  the 
downstream  droplet. 


SURFACE  TANGENTIAL  VELOCITY 


_  Lead  droplet,  t  s  1.0,  Re(t) »  92.2,  S.P.  at  141.4,  Local  max.  at  86.3 

_ Lead  droplet,  t  *  2.0,  R4(t)  =  85.9,  S.P.  at  148.9,  Local  max.  at  86.5 

_ Lead  droplet,  t  s  3.0,  Re(t)  =  80.7,  S.P.  at  161.2,  Local  max.  at  86.8 

_  Downstream  droplet,  t  =  1.0,  Re(t)  =  94.2,  S.P.  at  144.7,  Local  max.  at  87.2 

_  Downstream  droplet,  t  =  2.0,  Re(t)  =  90.3,  S.P.  at  157.6,  Local  max.  at  94.5 

_ Downstream  droplet,  t  =  3.0,  Re(t)  =  87.1,  S.P.  at  12.0,  Local  max.  at  104.7 


ANGULAR  POSITION  (DEGREE) 

Fuel :  N-DECANE 
Ambient  Temperature  >1000.  K 
Initial  Droplet  Temperature  ■  300.  K 
Initial  Reynolds  Number  ■  100. 

Initial R1  >  1j00  ,  R2>1j00  ,  D.8.0 

Figure  D- 14(c)  Surface  tangential  velocity  distribution  at  different  times  for  the  lead  droplet  and 
the  downstream  droplet. 


Lead  droplet,  ts  1.0,  Re(t)  =  92.2,  Local  min.  at  95.4 

Lead  droplet,  t  =  2.0,  R#(t)  =  85.9,  Local  min.  at  95.3 

Lead  droplet,  t  s  3.0,  Re(t)  =  80.7,  Local  min.  at  88.0 

Downstream  droplet,  ts  1.0,  Re(t)  =  94.2,  Local  min.  at  96.7 
Downstream  droplet,  t  *  2.0,  Re(t)  =  90.3,  Local  min.  at  97.3 
Downstream  dropleL  t  =  3.0,  Re(t)  =  87.1,  Local  min.  at  105.2 


0.  18.  36.  54.  72.  90.  108.  126.  144.  162.  180. 

I  ANGULAR  POSITION  (DEGREE) 

Fuel :  N-DECANE 

I  Ambient  Temperature  *1000.  K 

Initial  Droplet  Temperature  *  300.  K 
Initial  Reynolds  Number  *  100. 

|  Initial R1  ■  1.00  ,  R2 ■  1  j00  ,  D.8.0 

Figure  D- 14(d)  Surface  nondimension&l  pressure  distribution  at  different  times  for  the  lead  droplet 
|  and  the  downstream  droplet. 


- Load  droplet,  t  *  1.0,  Re(t)  *  92.2,  S.P.  at  141.4,  Local  min.  at  149.9 

- Lead  droplet,  t«  2.0,  Re(t)  =  85.9,  S.P.  at  148.9,  Local  min.  at  150.6 

- Lead  droplet,  t*  3.0,  Ra(t)  =  80.7,  S.P.  at  161.2,  Local  min.  at  151.0 

_ Downstream  droplet,  t  =  1.0,  Re(t)  =  94.2,  S.P.  at  144.7,  Local  min.  at  139.8 

_ Downstream  droplet,  t  =  2.0,  Re(t)  =  90.3,  S.P.  at  157.6,  Local  min.  at  0.0 

_ Downstream  droplet,  t  s  3.0,  Ra{t)s  87.1,  S.P.  at  12.0,  Local  min.  at  0.0 


ANGULAR  POSITION  (DEGREE) 

Fuel :  N-DECANE 
Ambient  Temperature  si 000.  K 
Initial  Droplet  Temperature  ■  300.  K 
Initial  Reynolds  Number  «  100. 

Initial  R1  ■  IjOO  ,  R2  ■  UOO  ,  D  >  84) 

Figure  D-14(e)  Local  Nusselt  number  distribution  at  different  times  for  the  lead  droplet  and  the 
downstream  droplet. 


LOCAL  SHERWOOD  NUMBER  Sh 


_  Lead  droplet,  t=  1.0,  Re(t)  =  02.2,  S.P.  at  141.4,  Local  min.  at  129.9 

_ Lead  droplet,  t  =  2.0,  Re(t)  =  85.9,  S.P.  at  148.9,  Local  min.  at  132.2 

_ Lead  droplet,  t  =  3.0,  Re(t)  =  80.7,  S.P.  at  161 J2,  Local  min.  at  132.1 

_  Downstream  droplet,  t  =  1.0,  Re(t)  =  94.2,  S.P.  at  144.7,  Local  min.  at  0.0 

_  Downstream  droplet,  t  a  2.0,  Re(t)  a  90.3,  S.P.  at  157.6,  Local  min.  at  18.0 

_ Downstream  droplet,  t  =  3.0,  Re(t)  =  87.1 ,  S.P.  at  12.0,  Local  min.  at  0.0 


ANGULAR  POSITION  (DEGREE) 

Fuel :  N-DECANE 
Ambient  Temperature  al 000.  K 
Initial  Droplet  Temperature  a  300.  K 
Initial  Reynolds  Number  a  100. 

Initial  R1  a  1.00  ,R2aU00  ,Di|j0 


Figure  D-  14(f)  Local  Sherwood  number  distribution  at  different  times  for  the  lead  droplet  and  the 
downstream  droplet. 


_  Lead  droplet,  Is  1.0,  Re(t)  =  92.2 

_ Lead droplet,  Is  2.0, Re(t)  =  85.9 

_ Lead  droplet,  ts  3.0,  Re(t)  =  80.7 

_  DownstrMm  droplet,  t  *  1.0,  Re(t)  =  94.2 

_ Downstream  droplet,  I  *  2.0,  Re(t)  =  90.3 

_ Downstream  droplet,  t  s  3.0,  Re(t)  =  87.1 


ANGULAR  POSITION  (DEGREE) 

Fuel :  N-DECANE 
Ambient  Temperature  s1000.  K 
Initial  Droplet  Temperature  *  300.  K 
Initial  Reynolds  Number  ■  100. 

Initial RIslJOO  ,  R2 * IjOO  ,  D«8 JO 


Figure  D-14(g)  Surface  temperature  distribution  at  different  times. 


DRAG  COEFFICIENTS 


_ 1.  Isolated  Droplet 

_ 2.  Lead  Droplet,  R1=1,  R2=1,  D=8,  Rel=100 

_ 3.  Downstream  Droplet,  R1=1,  R2=1,  D=8,  Re  1=1 00 


INSTANTANEOUS  REYNOLDS  NUMBER 


Figure  D-15  Time  variation  of  drag  coefficients  for  the  lead  droplet  and  the  downstream  droplet. 


DRAG  COMPONENTS 


_ 1.  Isolated  Droplet,  Pressure  Drag 

_ 2.  Lead  Droplet,  Pressure  Drag 

_ 3.  Downstream  Droplet,  Pressure  Drag 

_ 4.  Isolated  Droplet,  Friction  Drag 

_ 5.  Lead  Droplet,  Friction  Drag 

_ 6.  Downstream  Droplet,  Friction  Drag. 

_ 7.  Isolated  Droplet,  Thrust  Drag 

_ 6.  Lead  Droplet,  Thrust  Drag. 

_ 9.  Downstream  Droplet,  Thrust  Drag. 


GAS  HYDRODYNAMIC  DIFFUSION  TIME  SCALE 


Figure  D-16  Time  variation  of  three  drag  components  for  the  lead  droplet  and  the  downstream 
droplet. 


NUSSELT  NUMBER 


.  1.  Isolated  Droplet 

..  2.  Lead  Droplet,  R1=1,  R2=1,  D=8,  Rel=100 
_  3.  Downstream  Droplet,  R1=1,  R2=1,  D=8,  Re  1=1 00 
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GAS  HYDRODYNAMIC  DIFFUSION  TIME  SCALE 

Figure  D-17  Time  variation  of  average  Nusselt  numbers  for  the  lead  droplet  end  the  downstream 
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_ 1.  Isolated  Droplet 

_ 2.  Lead  Droplet,  R1=1,  R2=1,  D=8,  Rei=100 

_ 3.  Downstream  Droplet,  R1=1,  R2=1,  D=8,  Rel=100 


GAS  HYDRODYNAMIC  DIFFUSION  TIME  SCALE 


Figure  D-18  Time  variation  of  the  average  volumetric  temperature  for  the  lead  droplet  and  the 
downstream  droplet. 


SHERWOOD  NUMBER 


_ 1.  Isolated  Droplet 

_ 2.  Load  Droplet,  R1=1,  R2=1,  D=8,  Rel=100 

_ 3.  Downstream  Droplet,  R1=1,  R2=1,  D=8,  Rei=100 


GAS  HYDRODYNAMIC  DIFFUSION  TIME  SCALE 


Figure  D-19  Time  variation  of  average  Sherwood  numbers  for  the  lead  droplet  and  the  downstream 
droplet. 


C-C  SPACING  BETWEEN  TWO  DROPLETS 


_ 1.  Tg=1000  K,  H=300  K,  n-decane,  R1=1,  R2=1,  D=4 

_ 2.  Tg=1000  K,  Ti=300  K,  n-decane,  R1=1,  R2=1,  D=6 

_ 3.  Tg=1000  K,  Ti=300  K,  n-decane,  R1=1,  R2=1,  D=8 

_ 4.  Tg=1000  K,  Ti=300  K,  n-decane,  R1=1,  R2=1,  D=10 

_ 5.  Tg=1000  K,  Tl=300  K,  n-decane,  R1=1,  R2=1,  D=16 


GAS  HYDRODYNAMIC  DIFFUSION  TIME  SCALE 


Figure  D-20  Time  variation  of  droplet  (pacing  for  different  initial  droplet  (paring. 


DRAG  COEFFICIENTS 


_ 1.  Lead  Droplet,  R1=1,  R2=1,  D=4 

_ 2.  Lead  Droplet,  R1=1,  R2=1,  D=8 

_ 3.  Lead  Droplet,  R1=1,  R2=1,  D=16 

_  4.  Downstream  Droplet,  R1=1,  R2=1,  D=4 

_ 5.  Downstream  Droplet,  R1=1,  R2=1,  D=8 

_ 6.  Downstream  Droplet,  R1=1,  R2=1,D=16 
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INSTANTANEOUS  REYNOLDS  NUMBER 


Figure  D-21  Time  variation  of  drag  coefficients  for  different  initial  droplet  spacing. 


C-C  SPACING  BETWEEN  TWO  DROPLETS 


_  1.  Tg=1000  K,  Ti=300  K,  n-decane,  R1=1,  R2=0.25,  D=8 

_ 2.  Tg=1000  K,  Tl=300  K,  n-decane,  R1=1,  R  2=0.5,  D=8 

_ 3.  Tg=1000  K,  Ti=300  K,  n-decane,  R1=1,  R2=0.75,  D=8 

_ 4.  Tg=lOOO  K,  Tl=300  K,  n-decane,  R1=1,  R2=1,  D=8 

_ 5.  Tg=1000  K,  Tl=300  K,  n-decane,  R1=1,  R2=1.25,  D=8 


GAS  HYDRODYNAMIC  DIFFUSION  TIME  SCALE 


Figure  D-22  Time  variation  of  droplet  spacing*  for  different  downstream  droplet  radius. 


DRAG  COEFFICIENTS 


_ 1.  Usd  Droplet,  Riel,  R2=OJ25,  D=8 

_ 2.  Lead  Droplet,  R1=1,  R2=0.5,  D=8 

_ 3.  Lead  Droplet,  R1=1,  R2=0.75,  D=8 

_ 4.  Lead  Droplet,  R1=1,  R2=1n  D=8 

_ 5.  Downstream  Droplet,  R1=1,  R2=0.25,  D=8 

_ 6.  Downstream  Droplet,  R1=1,  R2=0.5,  D=8 

_ 7.  Downstream  Droplet,  Rial,  R2=0.75,  D=8 

_ 8.  Downstream  Droplet,  Rial,  R2=1,  D=8 


INSTANTANEOUS  REYNOLDS  NUMBER 


Figure  D-23  Time  variation  of  drag  coefficients  for  different  downstream  droplet  radius. 


C-C  SPACING  BETWEEN  TWO  DROPLETS 


_ 1.  Tg=1000  K,  Ti=300  K,  n-decane,  R1=1,  R2=1,  D=8,  Rel=50 

_ 2.  Tg=1000  K,  Tl=300  K,  n-decane,  R1=1,  R2=1,  D=8,  Re i =100 

_ 3.  Tgsiooo  K,  T1=300  K,  n-decane,  R1=1,  R2=1, 0=8,  Rel=125 


GAS  HYDRODYNAMIC  DIFFUSION  TIME  SCALE 


Figure  D-24  Time  variation  of  droplet  epadngs  for  different  initial  Reynolds  numbers. 


C-C  SPACING  BETWEEN  TWO  DROPLETS 


-  1.  Tg=1000  K,  Ti=300  K,  n-decane,  R1=1,  R2=0.7S,  D=8,  Re i =100 
..  2.  TgslOOO  K,  H=300  K,  n-decane,  R1=1,  R2=0.75,  D=8,  Rel=35 


GAS  HYDRODYNAMIC  DIFFUSION  TIME  SCALE 


Figure  D-25  Change  of  droplet-motion  pattern  from  coalescence  to  separation  by  changing  initial 
Reynolds  numbers. 


_ 1.  Tg=1000  K,  H=300  K,  n-decane,  R1=1,  R2=0.7,  D=4,  Re  1=100 

_ 2.  To =1000  K,  71=300  K,  n-decane,  R1=1,  R2=0.5,  D=4,  Rel=100 

_ 3.  Tb=1000  K,  Tl=300  K,  n-decane,  R1=1,  R2=0.7,  D=4,  Rei=50 

_ 4.  Tg=1000  K,  Ti=300  K,  n-decane,  Rial,  R2=0.5,  D=4,  Rei=50 


GAS  HYDRODYNAMIC  DIFFUSION  TIME  SCALE 


Figure  D-26  Time  variation  of  droplet  spadngs  for  different  downstream  droplet  radius  and  different 
initial  Reynolds  numbers. 


